Electromechanical model and simple numerical analysis for a piezoelectric vibration energy harvester considering nonlinear piezoelectricity,nonlinear damping,and self-powered synchronized switch circuit

Abstract

Establishment of an analytical technique that considers nonlinear piezoelectricity, nonlinear damping, and the connections to the self-powered synchronized switch circuit represents a challenging but practical issue in the development of piezoelectric vibration energy harvesters. The two-way coupled analysis method, which combines numerical software with a circuit simulator, can simulate the performance of the harvester, but it imposes a high computational load. We develop a numerical analysis technique that is more readily implemented to allow us to produce accurate predictions of the output power, displacement, and frequency response characteristics of this harvester with the switch circuit. First, we derive the electromechanical equation required for the miniature harvester, and we then determine the parameters in this equation using the harmonic balance method. Finally, the governing equation coupling the piezoelectric vibration energy harvester with the self-powered synchronized switch circuit is simplified, taking the form of a quartic equation with respect to the displacement. The simulation replicated a reduction in the displacement at higher resistive loads, wider variations in the resonance frequency, and the nonlinear stiffness softening characteristics that were observed experimentally. This paper provides a practical method for predicting the performance of the nonlinear piezoelectric vibration energy harvester with the synchronized switch circuit.

Keywords

Piezoelectric vibration energy harvester material nonlinearity synchronized switch harvesting on inductor self-powered circuit switch delay

1. Introduction

The piezoelectric vibration energy harvester (PVEH) is a promising power source for battery-less wireless sensor nodes because of its high power density and high adaptability to miniaturization. One major technical issue in the attempt to reduce the intervals between sensing data transmissions is maximization of the PVEH’s output power while maintaining miniaturization. To increase the efficiency of the PVEH’s circuit, approaches called synchronized switch harvesting on inductor (SSHI) or on a bias-flip circuit (Guyomar et al., 2005), which uses an active switch and an inductor, have been proposed. The parallel (Guyomar et al., 2005) and series (Badel et al., 2006) SSHI circuits flipped and magnified the piezoelectric voltage by switching to the inductor at voltages around the maximum/minimum piezoelectric voltage, thereby increasing the device’s harvesting power. The synchronous electric charge extraction (SECE) circuit (Lefeuvre et al., 2005) can reduce the load dependence of the output power by temporarily storing harvested energy in the inductor. Subsequently, the self-powered switching technique (Richard et al., 2007), which uses a proportion of the harvested energy to control the switch, was developed. Based on this technique, the self-powered SSHI (Lallart and Guyomar, 2008; Liang and Liao, 2012) and SECE (Zhao and Yang, 2015) were proposed. For analysis of these circuits, Lallart et al. (2012) and Chen et al. (2019) indicated the importance of switching delays induced by nonnegligible parasitic components in the discrete circuits. In addition to use of these self-powered switch circuits, use of a microcontroller unit with low energy consumption to control the switch represents another promising approach (Hara et al., 2019, 2020).

Unlike the simple waveform for a harvesting circuit consisting of a single resistive load, the piezoelectric voltage waveform for the switch circuits includes nonlinear (or precipitous) changes caused by switching, which makes it difficult to analyze the electromechanically coupled equation required. An analysis assuming the tip velocity of the PVEH and the external force are assumed to be in-phase (Badel et al., 2006; Guyomar et al., 2005; Lefeuvre et al., 2005) provides an effective methodology to determine the output power and the displacement at resonance, although the PVEH’s frequency response characteristics cannot be obtained because the assumption above does not hold outside the resonance. Shu et al. (2007) derived an analytical solution for the frequency response curves of a PVEH connected to SSHI circuits. They reported that the SSHI can only increase the output power significantly when the value of k_e²/ζ_m is much less than 1, where k_e and ζ_m are the effective electromechanical coupling coefficient and the mechanical damping ratio, respectively; otherwise, the SSHI cannot boost output power because the increased piezoelectric coupling force produced by the SSHI reduces the displacement and thus reduces the output power. To overcome this limitation of the SSHI, Makihara et al. (2015) and Makihara and Asahina (2017) proposed an intermittent switching technique to increase the output power by easing the increased piezoelectric coupling force. Morel et al. (2018) proposed the N-SECE technique with use of multiple and regenerative (or intermittent) switching to broaden the frequency bandwidth of a PVEH with a high electromechanical coupling coefficient. Similarly, the tunable (Lefeuvre et al., 2017) and frequency-tuning (Brenes et al., 2020a) SECE techniques can achieve wide frequency bands by using highly coupled states. A thorough review of the derived techniques can be found in the literature (Brenes et al., 2020b).

The electromechanically coupled equation based on linear piezoelectricity is widely accepted for use in analysis of the displacement and piezoelectric voltage characteristics of a PVEH connected to the aforementioned switch circuits. However, recent studies have shown that the material nonlinearities of the piezoelectric elements embedded in the structure should be included when modeling PVEHs, actuators, or piezoelectric attenuation devices to predict the frequency response characteristics accurately. Nonlinear piezoelectricity appears when high stress or strain is applied to piezoelectric elements. The nonlinear piezoelectricity primarily causes nonlinear stiffness softening, which affects the peak frequency, peak magnitude, and shape of the frequency response curve significantly, even when the displacement of the PVEH is small enough to ensure geometric linearity (Yang and Upadrashta, 2016). Through consideration of the nonlinearity in the relationship between the displacement and the excitation voltage of the piezoelectric actuator, Wagner and Hagedorn (2002) formulated the nonlinear enthalpy density and the nonlinear constitutive equation to derive the equation of motion for the actuator. The nonlinear terms of the enthalpy density were expressed as a polynomial form with respect to the strain. Mahmoodi et al. (2008) considered both the nonlinear piezoelectricity and the geometric nonlinearity when modeling the piezoelectrically actuated microcantilever sensor. By plotting a backbone curve that traced the peak values of the frequency response curves at different excitation levels, they showed that the degree of stiffness softening caused by the nonlinear piezoelectricity exceeded the degree of stiffness hardening caused by the geometric nonlinearity within the specific excitation level. Based on the polynomial model proposed by Wagner and Hagedorn (2002), an analysis of a PVEH connected to a single resistive load circuit (Meesala and Hajj, 2018; Stanton et al., 2010a, 2012), a finite element simulation of a simple unimorph PVEH connected to an SSHI circuit (Hegendörfer et al., 2022), and an analysis of a piezoelectric attenuation device connected to an autonomous SSD circuit (Asanuma and Komatsuzaki, 2020) were reported. In another viable model of nonlinear piezoelectricity, Goldschmidtboeing et al. (2011) considered the ferroelastic hysteresis phenomena observed within the strain-stress relationship of piezoelectric materials and derived the equation of motion by incorporating the hysteresis into the model. Leadenham and Erturk (2015) formulated the nonlinear enthalpy density based on the work of Goldschmidtboeing et al. (2011) and then derived a distributed parameter model for a piezoelectric actuator and a PVEH. Unlike the polynomial model proposed by Wagner and Hagedorn (2002), the nonlinear terms of the enthalpy density based on the ferroelastic hysteresis are dependent on the absolute value (or the signum function) of the strain, which is characterized by the hysteresis loop. Recent studies of a multimodal PVEH (Machado et al., 2020), a PVEH with segmented electrodes (Hu et al., 2021) connected to a single resistive circuit, and a PVEH connected to a diode bridge rectifier circuit (Leadenham and Erturk, 2020) validated the model proposed by Leadenham and Erturk (2015) by replicating the experimental results well. In addition to nonlinear stiffness softening, some studies also reported the importance of including nonlinear damping in the modeling process to explain the experimentally observed result that the displacement magnitude decreased when the excitation level increased. Several studies proposed models of nonlinear damping caused by a hysteresis damping (Goldschmidtboeing et al., 2011), nonconservative relations in nonlinear piezoelectricity (Stanton et al., 2012), or by structural energy dissipation (Asanuma and Komatsuzaki, 2020; Stanton et al., 2010b; Yang and Upadrashta, 2016). As methodologies to determine the nonlinear coefficients numerically, the harmonic balance method (Ito et al., 2020; Leadenham and Erturk, 2015; Stanton et al., 2010a; Wagner and Hagedorn, 2002), the method of multiple scales (Meesala and Hajj, 2018; Stanton et al., 2012), and the Bayesian updating technique (Poblete et al., 2021) were used. As noted by Shivashankar and Kandagal (2019) when determining the nonlinear coefficients, the applied nonlinear model should also be able to replicate the experimental results for the short-circuited frequency response curves of the displacement, in which the piezoelectric coupling force vanishes and thus only the nonlinear elastic terms affect the frequency response. In accordance with this observation, Asanuma (2021) reported that the model by Leadenham and Erturk (2015) provided a better fit to the experimental results for the frequency response curves of the piezoelectric attenuation device, while the model by Wagner and Hagedorn (2002) exhibited jumping phenomena caused by excessively strong nonlinear elastic terms that were not observed in the experiments.

The nonlinear piezoelectricity, the nonlinear damping behavior, the nonlinear change in the piezoelectric voltage caused by switching, the switching delay, and/or the phase difference between the displacement of the PVEH and the external force should be considered during modeling to provide accurate predictions of the output power, displacement, and frequency response characteristics of the PVEH when connected to the self-powered synchronized switch circuit. The two-way coupled analysis method, which combines numerical software with a circuit simulator (Asanuma et al., 2018, 2019), can simulate the real performance of a PVEH connected to synchronized switch circuits. However, use of the two-way coupled analysis involves a high computational load. In this study, we develop a more readily implemented numerical analysis technique that allows us to predict the output power, displacement, and frequency response characteristics of a PVEH when connected to a synchronized switch circuit. First, we derive the electromechanical equation for a miniature PVEH with a tip mass while considering both the nonlinear piezoelectricity and the structural discontinuity. Then, we determine the parameters for this equation using the harmonic balance method. Finally, we derive a governing equation to couple the harvester with the self-powered synchronized switch circuit and then verify that the solution to the governing equation replicates the experimental results.

2. PVEH structure and model

2.1. Modeling a miniature PVEH

We derive the electromechanical equation for the PVEH based on Hamilton’s principle and a Galerkin method. Figure 1 shows the structure of the miniature PVEH with a tip mass that is considered in this study. The PVEH consists of a beam (with a phosphor bronze strip for a spring), two tip masses (tough pitch copper), and two thin piezoelectric materials (lead zirconate titanate (PZT)-based soft ceramic C-93, Fuji Ceramics, Fujinomiya, Japan) that were bonded manually using epoxy adhesive on both sides. For simplicity, the beam and the tip masses were assumed to be made from the same copper alloy, thus meaning that they had the same mass density. The x, y, and z spatial axes were aligned with the length, width, and height directions of the PVEH, respectively. The terminals of the piezoelectric materials were connected in parallel. In the figure, b is the width of the PVEH, l_b is the length of the gap between the fixed edge of the PVEH and the left-side edge of the tip mass, l_t is half the length of the tip mass, l_p is the length of the piezoelectric material layer and the epoxy adhesive layer, and h_t, h_p, h_e, and h_b are the heights of the single tip mass, the piezoelectric material, the epoxy adhesive, and the beam, respectively. In a miniature harvester structure, the point mass approximation may not be applicable if the length of the tip mass is nonnegligible when compared with the beam length; therefore, in a manner similar to previous studies (Mak et al., 2011; Wang et al., 2020), this study considers a tip mass dimension to determine the kinetic energy. Additionally, this study incorporates an epoxy adhesive layer during modeling of the PVEH because the epoxy adhesive layer thickness increases the piezoelectric coupling coefficient when 2h_e is nonnegligible in comparison to h_b + h_p (Asanuma and Komatsuzaki, 2020). w_ext(x,t) and w(x,t) are the absolute displacement imposed externally by the shaker and the relative displacement from the stationary condition at position x on the beam, respectively. The material and structural parameters used in this work are summarized in Table 1. The material parameters of the custom-ordered piezoelectric material were evaluated based on the inspection report measurements provided by the manufacturer.

Figure 1.

(a) Structure of the miniature piezoelectric vibration energy harvester (PVEH) and (b) sectional view of the energy harvesting components.

Table 1.

Material and structural parameters.

$c_{11}^{b}$	110	GPa
$c_{11}^{p}$	63	GPa
$c_{11}^{e}$	2	GPa
e ₃₁	−15.2	C/m²
C _p	108	nF
ρ _c	8800	kg/m³
ρ _p	7910	kg/m³
ρ _e	1170	kg/m³
h _b	0.3	mm
h _p	0.2	mm
h _e	0.05	mm
h _t	4	mm
b	20	mm
l _b	33	mm
l _p	12	mm
l _t	2	mm

The kinetic energy T of the PVEH system can be obtained from the following equation.

\begin{matrix} T = \int_{0}^{l_{p}} \frac{1}{2} m_{eh} {(\overset{\cdot}{w} + {\overset{\cdot}{w}}_{ext})}^{2} dx + \\ \int_{l_{p}}^{l_{b}} [\frac{1}{2} m_{b} {(\overset{\cdot}{w} + {\overset{\cdot}{w}}_{ext})}^{2} + \frac{1}{2} M {(\overset{\cdot}{w} + {\overset{\cdot}{w}}_{ext})}^{2} δ (x - l_{b}) + M l_{t} \overset{\cdot}{w} \overset{\cdot}{w}' δ (x - l_{b}) + \frac{1}{2} J {\overset{\cdot}{w}'}^{2} δ (x - l_{b})] dx, \end{matrix}

(1)

where m_eh and m_b are mass per unit length of the energy harvesting components and the bare beam, respectively. M and J are the total tip mass and the mass’s inertia moment around the free end of the beam, respectively, δ is the Dirac delta function, and the dot and the prime above w and w_ext represent the time-derivative and the partial differential $\partial / \partial x$ , respectively. m_eh, m_b, and M are (ρ_ch_b+2 ρ _eh_e+2 ρ _ph_p)b, ρ_ch_bb, and 2ρ_c(2h_t+h_b)l_tb, where ρ_c, ρ_e, and ρ_p are the mass densities of the copper alloy, the epoxy adhesive, and the piezoelectric material, respectively. J was evaluated as follows:

J = \int_{l_{b}}^{l_{b} + 2 l_{t}} \int_{- h_{t} - h_{b} / 2}^{h_{t} + h_{b} / 2} ρ_{c} b [{(x - l_{b})}^{2} + z^{2}] dzdx .

(2)

The displacement magnitude of the PVEH used in this study is within the range of the geometric linearity, and we thus assume that the strain S₁ along the x-axis can be approximated using $- z w''$ . Based on this assumption, the potential energy U of the PVEH system is obtained from the following volume integration performed over the PVEH structure.

\begin{matrix} U = \int_{0}^{l_{p}} \int_{- b / 2}^{b / 2} \int_{- h_{b} / 2}^{h_{b} / 2} \frac{1}{2} c_{11}^{b} S_{1}^{2} dzdydx + \int_{l_{p}}^{l_{b}} \int_{- b / 2}^{b / 2} \int_{- h_{b} / 2}^{h_{b} / 2} \frac{1}{2} c_{11}^{b} S_{1}^{2} dzdydx \\ + \int_{0}^{l_{p}} \int_{- b / 2}^{b / 2} (\int_{h_{b} / 2}^{h_{b} / 2 + h_{a}} \frac{1}{2} c_{11}^{e} S_{1}^{2} dz + \int_{- h_{b} / 2 - h_{a}}^{- h_{b} / 2} \frac{1}{2} c_{11}^{e} S_{1}^{2} dz) dydx \\ + \int_{0}^{l_{p}} \int_{- b / 2}^{b / 2} (\int_{h_{b} / 2 + h_{a}}^{h_{b} / 2 + h_{a} + h_{p}} H_{p} dz + \int_{- h_{b} / 2 - h_{a} - h_{p}}^{- h_{b} / 2 - h_{a}} H_{p} dz) dydx, \end{matrix}

(3)

where $c_{11}^{b}$ and $c_{11}^{e}$ are the elastic stiffness constants of the beam and the epoxy adhesive, respectively. H_p is the electric enthalpy density based on the nonlinear piezoelectricity (Leadenham and Erturk, 2015) and is expressed as:

\begin{matrix} H_{p} = \frac{1}{2} c_{11}^{p} S_{1}^{2} + \frac{1}{3} c_{111}^{p} S_{1}^{2} | S_{1} | - e_{31} S_{1} E_{3} \\ - \frac{1}{2} e_{311} S_{1} | S_{1} | E_{3} - \frac{1}{2} ε_{33}^{S} {E_{3}}^{2}, \end{matrix}

(4)

where $c_{11}^{p}$ and $c_{111}^{p}$ are the linear and nonlinear elastic stiffness constants of the piezoelectric material, respectively, $e_{31}$ and $e_{311}$ are the linear and nonlinear piezoelectric constants, respectively, $ε_{33}^{S}$ is the permittivity at a constant strain, and E₃ is the electric field along the direction of the piezoelectric material’s thickness. The nonlinear parameters $c_{111}^{p}$ and $e_{311}$ are negative and positive, respectively, which means that the stiffness softening force and the electromechanical coupling both increase as S₁ increases. As will be explained in the subsection 2.2, these nonlinear parameters can be determined by fitting them to the experimental results. By substituting the approximated equation for S₁ into equation (3), the expression for U can be rewritten as:

\begin{matrix} U = \int_{l_{p}}^{l_{b}} \frac{1}{2} k_{0} {w''}^{2} dx \\ + \int_{0}^{l_{p}} [\frac{1}{2} k_{1} {w''}^{2} + \frac{1}{6} k_{2} {w''}^{2} | w'' | - \frac{1}{2} (2 θ_{1} w'' + θ_{2} w'' | w'' |) \overset{\cdot}{λ}] dx \\ - \frac{1}{2} C_{p} {\overset{\cdot}{λ}}^{2}, \end{matrix}

(5)

where k₀, k₁, k₂, θ₁, and θ₂ are all constant coefficients, λ(t) is the flux linkage, which is defined as a time-integral of the piezoelectric voltage v(t), and C_p is the capacitance of the parallel-connected piezoelectric materials. The five constant coefficients and C_p are written as follows:

k_{0} = \frac{1}{12} c_{11}^{b} {bh}_{b}^{3},

(6)

\begin{matrix} k_{1} = k_{0} + \frac{2}{3} c_{11}^{e} b [{(\frac{h_{b}}{2} + h_{e})}^{3} - {(\frac{h_{b}}{2})}^{3}] \\ + \frac{2}{3} c_{11}^{p} b [{(\frac{h_{b}}{2} + h_{e} + h_{p})}^{3} - {(\frac{h_{b}}{2} + h_{e})}^{3}], \end{matrix}

(7)

k_{2} = c_{111}^{p} b [{(\frac{h_{b}}{2} + h_{e} + h_{p})}^{4} - {(\frac{h_{b}}{2} + h_{e})}^{4}],

(8)

θ_{1} = e_{31} b (h_{b} + h_{p} + 2 h_{e}),

(9)

θ_{2} = \frac{1}{6} e_{311} b (3 h_{b}^{2} + 6 h_{p} h_{b} + 4 h_{p}^{2} + 12 h_{e}^{2} + 12 h_{e} h_{p} + 12 h_{e} h_{b}),

(10)

C_{p} = \frac{2 ε_{33}^{S} b l_{p}}{h_{p}} .

(11)

By considering the variables with respect to w and λ in equations (1) and (5), the variation in the Lagrangian L (=T−U) of the system can be expressed as (Leadenham and Erturk, 2015; Stanton et al., 2012):

δ L = (\frac{\partial L}{\partial \overset{\cdot}{w}}) δ \overset{\cdot}{w} + (\frac{\partial L}{\partial w''}) δ w'' + (\frac{\partial L}{\partial \overset{\cdot}{w}'}) δ \overset{\cdot}{w}' + (\frac{\partial L}{\partial \overset{\cdot}{λ}}) δ \overset{\cdot}{λ} .

(12)

In this study, the variation in the work by the nonconservative term is written as:

δ W_{nc} = - \int_{0}^{l_{b}} (d_{1} \overset{\cdot}{w} + d_{2} | \overset{\cdot}{w} | \overset{\cdot}{w}) dx δ w - i δ λ,

(13)

where d₁ and d₂ are the linear and nonlinear damping coefficients, respectively, and i is the electric current flowing to the circuit. A quadratic damping term with respect to the velocity is introduced in this study to replicate the previously reported experimental results in which the growth rate of the displacement decreased when the external acceleration increased (Asanuma and Komatsuzaki, 2020). The nonlinear damping can be attributed to an air-drag force acting on an elongated cantilever structure (Chen et al., 2017). Based on Hamilton’s principle (Leadenham and Erturk, 2015), δL and δW_nc satisfy the following relational equation between the given times t₁ and t₂.

\int_{t_{1}}^{t_{2}} (δ L + δ W_{nc}) dt = \int_{t_{1}}^{t_{2}} (δ T - δ U + δ W_{nc}) dt = 0 .

(14)

By solving equation (14) using the integration by parts approach, the electromechanical equations of the system can be derived and are expressed as follows:

\begin{matrix} [m_{b} + M δ (x - l_{b}) - M l_{t} δ' (x - l_{b})] \overset{\cdot\cdot}{w} + d_{1} \overset{\cdot}{w} + d_{2} | \overset{\cdot}{w} | \overset{\cdot}{w} + k_{0} w'''' \\ - J \overset{\cdot\cdot}{w}'' δ (x - l_{b}) - J \overset{\cdot\cdot}{w}' δ' (x - l_{b}) = - m_{b} {\overset{\cdot\cdot}{w}}_{ext} - M {\overset{\cdot\cdot}{w}}_{ext} δ (x - l_{b}) \\ {x | l_{p} < x < l_{b}}, \end{matrix}

(15)

\begin{matrix} m_{eh} \overset{\cdot\cdot}{w} + d_{1} \overset{\cdot}{w} + d_{2} | \overset{\cdot}{w} | \overset{\cdot}{w} + k_{2} (| w'' | w''')' - \\ {θ_{1} [δ' (x) - δ' (x - l_{p})] + θ_{2} w'''' sign (w'')} v = - m_{eh} {\overset{\cdot\cdot}{w}}_{ext} \\ {x | 0 < x < l_{p}}, \end{matrix}

(16)

\int_{0}^{l_{p}} [θ_{1} + θ_{2} | w'' |] \overset{\cdot}{w}'' dx + C_{p} \overset{\cdot}{v} + i = 0 {x | 0 < x < l_{p}},

(17)

where equations (15) and (16) represent the equations of motion for the bare beam (l_p < x < l_b) and the energy harvesting components (0 < x < l_p) of the PVEH, respectively, and equation (17) indicates the conservation law of electric current for the piezoelectric materials, where the sign is the signum function.

This study uses a discretized w(x,t) to simplify the electromechanical equation and considers a stiffness difference between the bare beam and the energy harvesting components. Therefore, w(x,t) is divided at the boundary between the bare beam and the energy harvesting components, and w(x,t) can be expanded as follows using a Galerkin method (Meesala and Hajj, 2018).

w (x, t) = \sum_{i = 1}^{N} ϕ_{b, i} (x) q_{i} (t) {x | l_{p} < x < l_{b}},

(18)

w (x, t) = \sum_{i = 1}^{N} ϕ_{eh, i} (x) q_{i} (t) {x | 0 < x < l_{p}},

(19)

where $ϕ_{b, i}$ and $ϕ_{eh, i}$ are the ith mode shapes of the bare beam and the energy harvesting components, respectively, q_i is the ith modal coordinate, and N is the number of modes. This study determines w(x,t) by solving a linear eigenvalue problem. Therefore, w(x,t) becomes the solution to the following eigenvalue equation.

m_{b} \overset{\cdot\cdot}{w} + k_{0} w'''' = 0 {x | l_{p} < x < l_{b}},

(20)

m_{eh} \overset{\cdot\cdot}{w} + k_{1} w'''' = 0 {x | 0 < x < l_{p}} .

(21)

By substituting equations (18) and (19) into equations (20) and (21), respectively, and by using the equation ${\overset{\cdot\cdot}{q}}_{i} = - ω_{i}^{2} q_{i}$ with respect to the ith angular eigenfrequency $ω_{i}$ , the eigenvalue equation can be re-expressed as:

ϕ''''_{b, i} = \frac{m_{b} ω_{i}^{2}}{k_{0}} ϕ_{b, i} = λ_{b, i}^{4} ϕ_{b, i} {x | l_{p} < x < l_{b}},

(22)

ϕ''''_{eh, i} = \frac{m_{eh} ω_{i}^{2}}{k_{1}} ϕ_{eh, i} = λ_{eh, i}^{4} ϕ_{eh, i} {x | 0 < x < l_{p}} .

(23)

As a general solution, $ϕ_{b, i}$ and $ϕ_{eh, i}$ are expressed in the following form.

\begin{matrix} ϕ_{b, i} (x) = A_{1} \cos (λ_{b, i} x) + A_{2} \sin (λ_{b, i} x) + A_{3} \cosh (λ_{b, i} x) \\ + A_{4} \sinh (λ_{b, i} x), \end{matrix}

(24)

\begin{matrix} ϕ_{eh, i} (x) = B_{1} \cos (λ_{eh, i} x) + B_{2} \sin (λ_{eh, i} x) + B_{3} \cosh (λ_{eh, i} x) \\ + B_{4} \sinh (λ_{eh, i} x), \end{matrix}

(25)

where A₁–A₄, B₁–B₄, $λ_{b, i}$ , and $λ_{eh, i}$ are constant coefficients. A₁–A₄, B₁–B₄, $λ_{b, i}$ , $λ_{eh, i}$ , and $ω_{i}$ were determined from the boundary conditions at the fixed end ( $ϕ_{eh, i} (0) = ϕ_{eh, i}' (0) = 0$ ), the continuity at x = l_p ( $ϕ_{eh, i} (l_{p}) = ϕ_{b, i} (l_{p})$ , $ϕ_{eh, i}' (l_{p}) = ϕ_{b, i}' (l_{p}),$ $k_{1} ϕ_{eh, i}'' (l_{p}) = k_{0} ϕ_{b, i}'' (l_{p})$ , $k_{1} ϕ'''_{eh, i} (l_{p}) = k_{0} ϕ'''_{b, i} (l_{p})$ ), and the balance of moment ( $k_{0} ϕ_{b, i}'' (l_{b}) = J ω_{i}^{2} ϕ_{b, i}' (l_{b}) + M l_{t} ω_{i}^{2} ϕ_{b, i} (l_{b})$ ) and the balance of force ( $k_{0} ϕ'''_{b, i} (l_{b}) = - M ω_{i}^{2} ϕ_{b, i} (l_{b}) - M l_{t} ω_{i}^{2} ϕ_{b, i}' (l_{b})$ ) at x = $l_{b}$ (Mak et al., 2011). Additionally, the following normalization was applied to determine the coefficients A₁–A₄ and B₁–B₄:

\begin{matrix} \int_{l_{p}}^{l_{b}} m_{b} ϕ_{b, i} ϕ_{b, j} dx + \int_{0}^{l_{p}} m_{eh} ϕ_{eh, i} ϕ_{eh, j} dx \\ + {(M ϕ_{b, i} ϕ_{b, j} + 2 M l_{t} ϕ_{b, i} ϕ_{b, j}' + J ϕ_{b, i}' ϕ_{b, j}') |}_{x = l_{b}} = δ_{ij}, \end{matrix}

(26)

where $δ_{ij}$ is the Kronecker delta.

Next, we combine the equations of motion given as equations (15) and (16) into a simpler single equation. The PVEH considered in this study is designed to generate power from the first vibrational mode. Hereafter, to simplify the notation, we omit the subscript mode number i from $ϕ_{b, i}$ , $ϕ_{eh, i}$ , and $q_{i}$ for the first vibrational mode. By multiplying the mode shapes $ϕ_{b}$ and $ϕ_{eh}$ for the first vibrational mode by both sides of equations (15) and (16), respectively, integrating the results with respect to x over a specified length, adding them together, and then simplifying using the normalization of equation (26) and the boundary conditions, the equation of motion with respect to q is derived and is expressed as follows:

\begin{matrix} \overset{\cdot\cdot}{q} + D_{1} \overset{\cdot}{q} + D_{2} | \overset{\cdot}{q} | \overset{\cdot}{q} + K_{1} q + K_{2} | q | q - (α_{1} + α_{2} | q |) v \\ = - M_{f} {\overset{\cdot\cdot}{w}}_{ext} . \end{matrix}

(27)

Similarly, by considering the first vibrational mode only, equation (17) can be re-expressed as follows:

(α_{1} + α_{2} | q |) \overset{\cdot}{q} + C_{p} \overset{\cdot}{v} + i = 0 .

(28)

The parameters above are expressed as follows:

D_{1} = \int_{l_{p}}^{l_{b}} d_{1} ϕ_{b}^{2} dx + \int_{0}^{l_{p}} d_{1} ϕ_{eh}^{2} dx,

(29)

D_{2} = \int_{l_{p}}^{l_{b}} d_{2} ϕ_{b}^{3} dx + \int_{0}^{l_{p}} d_{2} ϕ_{eh}^{3} dx,

(30)

K_{1} = ω_{1}^{2} = \int_{l_{p}}^{l_{b}} k_{0} {ϕ_{b}''}^{2} dx + \int_{0}^{l_{p}} k_{1} {ϕ_{eh}''}^{2} dx,

(31)

K_{2} = \int_{0}^{l_{p}} \frac{1}{2} k_{2} | ϕ_{eh}'' | {ϕ_{eh}''}^{2} dx,

(32)

α_{1} = {θ_{1} ϕ_{eh}' |}_{x = l_{p}},

(33)

α_{2} = \int_{0}^{l_{p}} θ_{2} | ϕ_{eh}'' | ϕ_{eh}'' dx,

(34)

M_{f} = \int_{l_{p}}^{l_{b}} m_{b} ϕ_{b} dx + \int_{0}^{l_{p}} m_{eh} ϕ_{eh} dx + {M ϕ_{b} |}_{x = l_{b}} .

(35)

2.2. Determination of parameters

Based on the approaches of previous studies (Ito et al., 2020; Leadenham and Erturk, 2015; Shivashankar and Kandagal, 2019), the harmonic balance method (HBM) was used to determine the parameters of equations (27) and (28) systematically. Following the important observations reported in previous studies (Asanuma, 2021; Shivashankar and Kandagal, 2019), the peak values of the frequency response characteristics for both short and open circuits were used to determine the parameters. q(t) and v(t) are assumed to be expressed using the following trigonometric functions under the harmonic excitation ${\overset{\cdot\cdot}{w}}_{ext} = A \cos ω t$ :

q (t) = q_{0} \cos (ω t + μ) = q_{1} \cos ω t - q_{2} \sin ω t,

(36)

v (t) = v_{0} \cos (ω t + ν) = v_{1} \cos ω t - v_{2} \sin ω t,

(37)

where ω is the angular frequency of the external excitation, q₀–q₂ and v₀–v₂ are the amplitudes, which have relationships based on the addition theorem of the trigonometric function (q₁ = q₀cosμ, q₂ = q₀sinμ, v₁ = v₀cosν, and v₂ = v₀sinν), and μ and ν are the corresponding phase differences with respect to ${\overset{\cdot\cdot}{w}}_{ext}$ , respectively. By substituting q(t) and v(t) into equations (27) and (28), considering the open-circuit condition (i.e. the case where i = 0 in equation (28)), and applying the first-order approximations of the harmonic terms, the following simultaneous equations for the open circuit case are derived.

\begin{matrix} - q_{1} ω^{2} - D_{1} q_{2} ω - \frac{8}{3 π} D_{2} q_{0} q_{2} ω^{2} + K_{1} q_{1} + \frac{8}{3 π} K_{2} q_{0} q_{1} \\ - α_{1} v_{1} - \frac{4 α_{2}}{3 π q_{0}} [(q_{0}^{2} + q_{1}^{2}) v_{1} + q_{1} q_{2} v_{2}] + M_{f} A = 0, \end{matrix}

(38)

\begin{matrix} q_{2} ω^{2} - D_{1} q_{1} ω - \frac{8}{3 π} D_{2} q_{0} q_{1} ω^{2} - K_{1} q_{2} - \frac{8}{3 π} K_{2} q_{0} q_{2} \\ + α_{1} v_{2} + \frac{4 α_{2}}{3 π q_{0}} [(q_{0}^{2} + q_{2}^{2}) v_{2} + q_{1} q_{2} v_{1}] = 0, \end{matrix}

(39)

- C_{p} v_{2} ω - α_{1} q_{2} ω - \frac{4}{3 π} α_{2} q_{0} q_{2} ω = 0,

(40)

- C_{p} v_{1} ω - α_{1} q_{1} ω - \frac{4}{3 π} α_{2} q_{0} q_{1} ω = 0 .

(41)

The simultaneous equations for the short circuit case can be derived by setting v₁ and v₂ to zero in equations (38) and (39). The parameters D₁, D₂, K₂, and α₂ were obtained numerically by fitting the simultaneous equation to the peak values of the experimentally-obtained frequency response curves. We used MATLAB^® (MathWorks) and its nonlinear least-squares solver to perform the fitting analysis (Asanuma, 2021). The parameters obtained from this analysis are listed in Table 2. Figure 2 shows the simulated (solid lines) and experimental (dots) results for the typical frequency response curves for the displacement and the open-circuit voltage v_oc at four different acceleration magnitudes (0.5, 1, 3, and 5 m/s²), along with the resonance frequency f_r, the displacement, and the open-circuit voltage at resonance versus the acceleration magnitudes. Increments in the frequency of 0.1 Hz were used in the frequency response curves around the resonance. Note that, as shown in Figure 2(e), the simulation results obtained for the short and open circuit cases were indistinguishable because the displacement magnitude at resonance was identical in the two cases. Overall, as shown in Figure 2, the simulation replicated the nonlinear stiffness softening behavior, the location of the resonance frequencies, and the magnitudes of the displacements and the voltages obtained in the experiments.

Table 2.

Parameters obtained from the fitting analysis.

M _f	7.1 × 10⁻²	kg^1/2
D ₁	2.86	s⁻¹
D ₂	69.2	kg^−1/2 m⁻¹
K ₁	1.59 × 10⁵	s⁻²
K ₂	−5.97 × 10⁷	kg^−1/2 m⁻¹ s⁻²
α ₁	−1.53 × 10⁻²	C kg^−1/2 m⁻¹
α ₂	−6.86 × 10⁻²	C kg⁻¹ m⁻²

Figure 2.

Simulated (solid lines) and experimental (dots) results for the typical frequency response curves of the displacements for the (a) open (OC) and (b) short (SC) circuits and for (c) the open-circuit voltage v_oc at four different acceleration magnitudes (0.5, 1, 2, 3, and 5 m/s²), and (d) the resonance frequency f_r, (e) the displacement, and (f) the open-circuit voltage at f_r versus the acceleration magnitudes.

2.3. Analysis coupled with circuits

The magnitude of q(t) was at most of the order of 10⁻⁴ within the excitation level range used in this study, and the nonlinear electromechanical coupling of the designed harvester, where the piezoelectric materials did not cover the beam fully, was low. Therefore, the magnitude of the term α₂|q| in equations (27) and (28) was negligibly small when compared with the magnitude of α₁. Accordingly, these two equations were simplified as follows:

\overset{\cdot\cdot}{q} + D_{1} \overset{\cdot}{q} + D_{2} | \overset{\cdot}{q} | \overset{\cdot}{q} + K_{1} q + K_{2} | q | q - α_{1} v = - M_{f} {\overset{\cdot\cdot}{w}}_{ext} .

(42)

α_{1} \overset{\cdot}{q} + C_{p} \overset{\cdot}{v} + i = 0 .

(43)

Figure 3 shows the standard (STD) circuit and the self-powered parallel SSHI (SP-SSHI) circuit (Richard et al., 2007) considered in this work. The STD and SP-SSHI were selected as classical but widely accepted circuits in the community working on the PVEHs. In addition, within this study using a low-coupled PVEH, the SP-SSHI is the best candidate in terms of power increase because the SECE or series SSHI may create a larger piezoelectric coupling (or damping) force and thus decrease the output power (Zhao and Yang, 2017). The STD circuit consists of a full bridge diode rectifier (BAT 48, STMicroelectronics), a smoothing capacitor C_s (10 μF), and a resistive load R_L to evaluate the output power. v_DC is the voltage across the R_L-C_s parallel circuit. The SP-SSHI circuit uses a self-powered switch circuit and an inductor L_i (47 mH) with an intrinsic resistor r, which were connected in parallel to the components that were used in the STD circuit. The detailed workflow for the self-powered switch can be found in the literature (Liang and Liao, 2012), and thus the lengthy explanation of this workflow is omitted here for brevity. The self-powered switch circuit was built using pnp (Tr₁ and Tr₄; NSS1C200L, onsemi) and npn (Tr₂ and Tr₃; MMBT8099LT1G, onsemi) bipolar transistors, Schottky diodes (BAT 48, STMicroelectronics), base resistors denoted by r_b (3.3 kΩ), and envelope-detecting sub-circuits consisting of a resistor r_ed (300 kΩ) and a capacitor C_ed (1 nF). When compared with the corresponding integrated circuits (Li et al., 2022), the self-powered switch constructed from discrete circuits can be designed easily and is resistant to high voltages/currents. However, the switching delay, which is caused by parasitic resistors and capacitors among the discrete circuit components, affects the harvester’s performance and increases the complexity of the coupled analysis.

Figure 3.

(a) Standard (STD) circuit and (b) self-powered parallel SSHI (SP-SSHI) circuit (Richard et al., 2007).

Next, we consider the coupled analysis between the nonlinear model of the PVEH and the STD and SP-SSHI circuits based on use of the analytical technique developed for the linear model of the PVEH (Shu et al., 2007; Zhang et al., 2021). Figure 4 shows schematics of typical waveforms for the modal coordinate q(t) and the piezoelectric voltage v(t) for (a) the STD circuit and (b) the SP-SSHI circuit. This study considers a specific waveform for the switch circuit to realize more accurate prediction of the frequency response characteristics when compared with the results obtained from previous studies, in which the waveform was approximated by using the signum function (Asanuma, 2021; Asanuma and Yamauchi, 2022).

Figure 4.

Schematics of waveforms for the modal coordinate q(t) and the piezoelectric voltage v(t) for (a) the STD circuit and (b) the SP-SSHI circuit.

First, we analyze the frequency response characteristics of the STD circuit by referring to the waveforms shown in Figure 4(a), and we then explain the corresponding analysis of the SP-SSHI circuit characteristics in a similar manner. When an external excitation, which is expressed as a sinusoidal function here, is applied to the PVEH, the terms $- M_{f} {\overset{\cdot\cdot}{w}}_{ext} (t)$ and q(t) in equation (42) can also be expressed as sinusoidal functions:

- M_{f} {\overset{\cdot\cdot}{w}}_{ext} = F (t) = F_{0} \sin ω t,

(44)

q (t) = q_{0} \sin (ω t - τ),

(45)

where F₀ and q₀ are the amplitudes, ω is the angular frequency of the external excitation, and τ is the phase difference with respect to the external excitation. As shown in Figure 4(a), the time interval between t₁ and t₃ is equal to one-half of the vibrational period T of q, during which q varies from a minimum (q(t₁) = −q₀) to a maximum (q(t₃) = q₀). The electric current outside the piezoelectric materials remains at zero when |v(t)| is less than v_c. When the piezoelectric voltage exceeds v_c at time t₂, the electric current begins to flow through the diode rectifier into the R_L-C_s parallel circuit. In reality, v_c is equal to the sum of v_DC and 2v_d, where v_d is the voltage drop across a single diode. To simplify the numerical analysis, it is assumed in this study that v_DC is much larger than 2v_d and is thus equal to v_c. Additionally, when the capacitance of the smoothing capacitor is large enough to allow the ripple of v_DC to be ignored and the time constant of R_LC_s is much greater than the vibrational period of q, v_DC is assumed to remain constant. Based on these assumptions, the electrical current i outside the piezoelectric materials satisfies the following equation (Shu et al., 2007):

\int_{t_{1}}^{t_{3}} idt = \frac{T v_{c}}{2 R_{L}} .

(46)

By solving equation (46) with reference to the waveforms shown in Figure 4(a), the following relational equation is derived with respect to v_c:

v_{c} = - \frac{2 α_{1} R_{L} ω}{π + 2 R_{L} C_{p} ω} q_{0},

(47)

where T/2 was rewritten as π/ω. Subsequently, we consider the energy balance of the system over the half vibration period from t₁ to t₃. By integrating over the half vibrational period after multiplying both sides of equation (42) by $\overset{\cdot}{q}$ , the following integral equation was derived:

\begin{matrix} \int_{t_{1}}^{t_{3}} \overset{\cdot\cdot}{q} \overset{\cdot}{q} dt + D_{1} \int_{t_{1}}^{t_{3}} {\overset{\cdot}{q}}^{2} dt + D_{2} \int_{t_{1}}^{t_{3}} | \overset{\cdot}{q} | {\overset{\cdot}{q}}^{2} dt + K_{1} \int_{t_{1}}^{t_{3}} q \overset{\cdot}{q} dt \\ + K_{2} \int_{t_{1}}^{t_{3}} | q | q \overset{\cdot}{q} dt - α_{1} \int_{t_{1}}^{t_{3}} v \overset{\cdot}{q} dt = \int_{t_{1}}^{t_{3}} F \overset{\cdot}{q} dt . \end{matrix}

(48)

By referring to the waveforms shown in Figure 4(a) and substituting in equations (43)–(45), equation (48) can be calculated as follows:

\frac{8}{3 π} D_{2} ω^{2} q_{0}^{2} + D_{1} ω q_{0} + \frac{8 α_{1}^{2} R_{L} ω}{{(π + 2 R_{L} C_{p} ω)}^{2}} q_{0} = F_{0} \sin τ .

(49)

We must now eliminate the phase difference τ on the right side of equation (49) to derive an equation in which the only independent variable is q₀. To derive another equation with respect to the phase difference τ (Shu et al., 2007), we calculated the time-derivative of equation (42), followed by calculation of the time-integral over the half vibrational period. The resulting equation is expressed as:

\begin{array}{l} \int_{t_{1}}^{t_{3}} \overset{⃛}{q} d t + D_{1} \int_{t_{1}}^{t_{3}} \overset{\cdot\cdot}{q} d t + 2 D_{2} \int_{t_{1}}^{t_{3}} | \dot{q} | \overset{\cdot\cdot}{q} d t + K_{1} \int_{t_{1}}^{t_{3}} \dot{q} d t \\ + 2 K_{2} \int_{t_{1}}^{t_{3}} | q | \dot{q} d t - α_{1} \int_{t_{1}}^{t_{3}} \dot{v} d t = \int_{t_{1}}^{t_{3}} \dot{F} d t . \end{array}

(50)

Similarly, by referring to the waveforms shown in Figure 4(a) and substituting in equations (43)–(45), equation (50) can be calculated as follows:

K_{2} q_{0}^{2} + K_{1} q_{0} - ω^{2} q_{0} + \frac{2 α_{1}^{2} R_{L} ω}{π + 2 R_{L} C_{p} ω} q_{0} = F_{0} \cos τ .

(51)

Therefore, by adding the squares of equations (49) and (51) and eliminating τ, the following equation, which is dependent on q₀ alone, was derived:

\begin{matrix} {(\frac{8}{3 π} D_{2} ω^{2} q_{0}^{2} + D_{1} ω q_{0} + \frac{8 α_{1}^{2} R_{L} ω}{{(π + 2 R_{L} C_{p} ω)}^{2}} q_{0})}^{2} \\ + {(K_{2} q_{0}^{2} + K_{1} q_{0} - ω^{2} q_{0} + \frac{2 α_{1}^{2} R_{L} ω}{π + 2 R_{L} C_{p} ω} q_{0})}^{2} = F_{0}^{2} . \end{matrix}

(52)

Second, we derive a similar q₀-dependent equation for the SP-SSHI circuit while considering the switching delay φ (Lallart et al., 2012). As shown in Figure 4(b), q(t) reaches a minimum (−q0) at t₀ and a maximum (q₀) at t₃. In this analysis, the time interval from t₁ to t₄ was defined as the half vibrational period T/2 because of the discontinuity of v(t) caused by switching. In addition, q(t₁) and q(t₄) are expressed as −q₀cosφ and q₀cosφ, respectively. The electric current i can flow through the diode rectifier into the R_L-C_S parallel circuit during the period from t₂ to t₃. v_M and v_m are the magnitudes of the piezoelectric voltage before and after switching, and the ratio of v_m/v_M is called the voltage inversion coefficient γ. By considering the voltage difference between v(t₃) and v(t₄), as shown in Figure 4(b), v_M can be expressed as follows:

v_{M} = v_{c} - [- \frac{α_{1}}{C_{p}} (q_{0} - q_{0} \cos ϕ)] = v_{c} + \frac{α_{1} q_{0}}{C_{p}} (1 - \cos ϕ) .

(53)

Similarly, by calculating equation (46) over the half period from t₁ to t₄ by substitution of equation (53), v_c can be expressed as follows:

v_{c} = A (ω) q_{0} = - \frac{α_{1} R_{L} [1 + \cos ϕ - γ (1 - \cos ϕ)] ω}{π + R_{L} C_{p} (1 - γ) ω} q_{0},

(54)

where A(ω) is the coefficient that is dependent on ω. Using a similar approach to the methodologies applied when performing the analysis for the STD circuit, by calculating and simplifying integral equations (48) and (50) over the half period from t₁ to t₄, the following two equations can be derived.

\begin{matrix} \frac{8}{3 π} D_{2} ω^{2} q_{0}^{2} + D_{1} ω q_{0} + \frac{2}{R_{L} ω} A^{2} q_{0} + \frac{C_{p}}{π} (1 - γ^{2}) \\ {[A + \frac{α_{1}}{C_{p}} (1 - \cos ϕ)]}^{2} q_{0} = B_{1} (ω, q_{0}) = F_{0} \sin τ, \end{matrix}

(55)

\begin{matrix} (K_{2} co s^{2} ϕ - D_{2} ω^{2} si n^{2} ϕ) q_{0}^{2} \\ + (K_{1} \cos ϕ + \frac{α_{1}^{2}}{C_{p}} \cos ϕ - M ω^{2} \cos ϕ - D_{1} ω \sin ϕ + \frac{π α_{1} A}{2 R_{L} C_{p} ω}) q_{0} \\ = B_{2} (ω, q_{0}) = F_{0} \cos (τ + ϕ), \end{matrix}

(56)

where B₁(ω, q₀) and B₂(ω, q₀), which represent the left-hand sides of equations (55) and (56), respectively, are the functions that are dependent on ω and q₀. By adding the squares of equations (55) and (56) to eliminate τ, the following equation was derived:

B_{1}^{2} + {(B_{1} \tan ϕ + \frac{B_{2}}{\cos ϕ})}^{2} = F_{0}^{2} .

(57)

q ₀ can be determined for the STD and SP-SSHI circuits by solving equations (52) and (57), respectively, after substituting in specific values for the parameters. v_c was then calculated for the STD and SP-SSHI circuits from equations (47) and (54), respectively, by substitution of the resulting q₀. The output power P_out was calculated as v_c²/R_L. Note that, by expanding and simplifying, equations (52) and (57) finally form the required quartic equation with respect to q₀, which can be solved analytically using Ferrari’s method (Weisstein, 2022). The analytical solution has quite a lengthy expression; therefore, we determined q₀ numerically in this study, based on considerations that the value of q₀ was real, positive, and close to the value determined at one frequency step ahead.

The tip displacement of the PVEH must be small enough to satisfy the condition of geometric linearity because we derived the quartic equation based on the geometric linearity. Figure 5 shows the results of a finite element analysis (FEA) of the restoring force versus the displacement of the PVEH based on the geometric nonlinearity and linearity. Comsol Multiphysics^® FEA software was used to perform this analysis. The structural model used in the FEA was identical to the structure shown in Figure 1. The displacement at the boundary between the mass and the beam was evaluated by applying a static force. The FEA indicated that the two results showed no differences at displacements smaller than 2 mm, whereas the nonlinear increase in stiffness caused by the geometric nonlinearity was clearly observed at displacements greater than 4 mm. Therefore, geometric linearity was ensured at displacements of less than 2 mm, and we thus proceeded with the numerical analysis while considering a geometrically linear range.

Figure 5.

Results of FEA of the restoring force versus displacement of the PVEH based on the geometric nonlinearity and the geometric linearity.

Figure 6 shows the numerical simulation results for the output power P_out, the displacement, and the resonance frequency f_r versus the resistive load R_L and versus the acceleration applied externally for the (a), (c), (e) STD and (b), (d), (f) SP-SSHI circuits. The output power and displacement values were acquired at the resonance frequency. The open circles connected by a white solid line represent the points at which the maximum of the output power was obtained for each given acceleration. The displacement was evaluated at the boundary between the mass and the beam, that is, where x = l_b. For simplicity and consistency, the inversion coefficient γ (=0.66) and the switching delay φ (=0.42 rad) were set to be constant in the analysis; these values are empirically attainable at high acceleration. As shown in Figure 6(a) and (c), the maximum output power for the STD circuit was obtained at the minimum displacement, and R_L remained approximately constant at the maximum output power. In contrast, as shown in Figure 6(b) and (d), the displacement of the SP-SSHI circuit continued to decrease as R_L increased, and R_L increased slightly at the maximum output power from 20 to 40 kΩ as the acceleration increased. This occurred because the SP-SSHI’s piezoelectric coupling force (−α₁v in equation (42)), which increases at larger values of R_L, suppressed the displacement significantly and thus reduced the output power, particularly at lower accelerations. Although the simulation located the maximum output power value successfully, we assume that experimental identification of a true maximum point may be difficult because, as illustrated in Figure 6(a) and (b), the shape of the distribution near the peak was rounded at R_L values between 20 and 60 kΩ. As shown in Figure 6(e) and (f), the resonance frequency for the SP-SSHI varied more widely when compared with that of the STD circuit, because the reduced displacement at larger R_L values reduced the stiffness softening caused by the nonlinear piezoelectricity, and connection to the SP-SSHI increased the stiffness, as shown in equation (56).

Figure 6.

Numerical simulations of output power P_out, displacement, and resonance frequency f_r versus the resistive load R_L and versus the externally applied acceleration for the (a, c, and e) STD circuit and the (b, d, and f) SP-SSHI circuit. The open circles connected using a white solid line represent the points at which the maximum output power was obtained at a given acceleration.

Figure 7 shows the results of numerical simulation of the output power P_out and the displacement versus the frequency and versus the acceleration that was applied externally for the (a), (c) STD circuit and (b), (d) the SP-SSHI circuit. To enable a fair comparison, the same resistive load of 36 kΩ was used in this case. Similar to Figure 6, the open circles connected using a white solid line represent the maximum output power values at each given acceleration. As shown in Figure 7, the resonance frequency shifted toward a lower frequency when the acceleration increased, and its asymmetric distribution with respect to the peak was clearly observed. Therefore, the stiffness softening behavior caused by the nonlinear piezoelectricity affects the frequency response characteristics significantly, and thus this softening cannot be neglected in the analysis. When compared with the STD circuit, use of the SP-SSHI circuit provided an increase in the output power, a reduction in the displacement, and alleviated the lower shift in the resonance frequency.

Figure 7.

Numerical simulation results for output power P_out and displacement versus frequency and versus the acceleration applied externally for (a and c) the STD circuit and (b and d) the SSHI circuit. The open circles connected using a white solid line represent the points at which the maximum output powers at a given acceleration were obtained.

3. Experimental results

In this section, we validate the established numerical analysis method via comparison with the experimental results. Figure 8 shows an experimental setup for measurement of the voltages and displacements of the piezoelectric harvester when connected to the STD and SP-SSHI circuits. A shaker system (shaker: SW-2015; controller: APD-200FCG; Asahi Seisakusyo) was used to apply a sinusoidal vibration, which was oriented perpendicular to the longitudinal direction of the piezoelectric harvester. Similar to the approach used in the simulation, the displacement was evaluated close to the boundary between the mass and the beam using a laser displacement sensor system (LK-H150 and LK-G5000, Keyence). The signals for the piezoelectric voltage v, the DC voltage v_DC across the resistive load, the displacement, and the acceleration applied by the shaker were collected through a four-channel differential oscilloscope (HS6-DIFF USB oscilloscope, TiePie Engineering), and the signals were then processed using a personal computer (PC). As shown in Figure 9, the inversion coefficient γ and the switching delay φ exhibited dependences on the acceleration (or the piezoelectric voltage) during the experiments because of the voltage-dependent parasitic elements present in the SP-SSHI circuit (Asanuma and Komatsuzaki, 2020). The parameter values of γ and φ obtained in this step were used in the numerical analysis. The other parameters used in the numerical analysis were found in Tables 1 and 2. All parameters of the PVEH are the same value as used in the numerical analysis.

Figure 8.

Experimental set-up for measurement of the voltages and displacement of the piezo harvester.

Figure 9.

Experimental results for the acceleration dependences of (a) the inversion coefficient γ and (b) the switching delay φ.

Figure 10 shows the resistive load R_L dependences of the output power P_out, the displacement, and the resonance frequency f_r at three different acceleration magnitudes (i.e. 1, 3, and 5 m/s²). Figure 11 shows the acceleration dependences of the output power P_out, the displacement, and the resonance frequency f_r. Figure 12 shows typical frequency response curves acquired at the three different acceleration magnitudes (1, 3, and 5 m/s²). To ensure a fair comparison, the same resistive load of 36 kΩ was used during evaluation of both the acceleration dependence and the frequency response. As predicted in subsection 2.3, the reduction in the displacement at larger values of R_L, a wider variation in the resonance frequency, and nonlinear stiffness softening were all observed experimentally, and the simulations replicated these phenomena well. Therefore, by using the established numerical analysis process, we can predict the electromechanical properties of a piezoelectric harvester connected to the SP-SSHI circuit with a lower computational load when compared with the two-way coupled analysis method using a combination of numerical software and a circuit simulator. For the SP-SSHI, a discrepancy between the experimental and numerical results was found at a larger R_L and/or a high acceleration of 5 m/s². The discrepancy for the SP-SSHI may be attributed to a higher vibrational mode induced by an increased piezoelectric coupling (or damping) force. After the switch closes at a large piezoelectric voltage (or, a high acceleration and/or a large R_L), a high current flows to the circuit and thus the impulsive damping force acts on the PVEH. The impulsive damping force may induce a higher vibrational mode, which triggers an unwanted switching and a reverse current. Figure 13 shows experimental results of waveforms of the piezoelectric voltage v, the DC voltage v_DC across the resistive load R_L, and the displacement for (a) R_L of 36 kΩ and (b) R_L of 1 MΩ at the acceleration of 5 m/s². As shown in Figure 13, a higher vibrational mode affected the displacement waveforms, and the unwanted weak switching was clearly observed as a ripple in the piezoelectric voltage waveform at R_L of 1 MΩ. Therefore, the presence of the higher vibrational mode may result in the discrepancy.

Figure 10.

Simulated and experimental results for the resistive load R_L dependences of (a, d and g) the output power P_out, (b, e and h) the displacement, and (c, f and i) the resonance frequency f_r at three different acceleration magnitudes (1, 3, and 5 m/s²).

Figure 11.

Simulated and experimental results for (a) the acceleration dependence of the output power P_out, (b) the displacement, and (c) the resonance frequency f_r at a resistive load of 36 kΩ.

Figure 12.

Simulated (solid lines) and experimental (dots) results for typical frequency response curves of (a and b) the output power P_out and (c and d) the displacement at three different acceleration magnitudes (1, 3, and 5 m/s²).

Figure 13.

Experimental results of waveforms of the piezoelectric voltage v, the DC voltage v_DC across the resistive load R_L, and the displacement for (a) R_L of 36 kΩ and (b) R_L of 1 MΩ at the acceleration of 5 m/s².

For simplicity, the voltage drops in the diode rectifier (Leadenham and Erturk, 2020; Yang et al., 2020) and the dielectric and piezoelectric losses of the piezoelectric materials (Gibus et al., 2022) which become significant for a PVEH with a high electromechanical coupling at a large piezoelectric voltage were not considered in this study. For accurate predictions, these factors should be incorporated into the model depending on the electromechanical coupling states, the magnitude of piezoelectric voltage, and/or acceleration level.

4. Conclusion

To replace the existing two-way coupled analysis method, which involves a high computational load, we have developed a numerical analysis technique that is more easily implemented to allow us to predict output power, displacement, and frequency response characteristics for a PVEH connected to a self-powered synchronized switch circuit. First, we derived the electromechanical equation required for a miniature PVEH, and we then determined the parameters in this equation using the harmonic balance method. Finally, the governing equation coupling the PVEH with the self-powered synchronized switch circuit can be simplified using the form of a quartic equation with respect to the displacement. The simulations successfully replicated the experimentally-observed phenomena, including a reduction in the displacement at higher resistive loads, a wider variation of the resonance frequency value, and nonlinear stiffness softening. Therefore, using the established numerical analysis technique, we can predict the electromechanical properties of a PVEH connected to a self-powered synchronized switch circuit with a reduced computational load. This study revealed that the piezoelectric coupling force increased by the self-powered synchronized switch circuit alleviates the nonlinear stiffness softening caused by the nonlinear piezoelectricity, thus those two factors cannot be ignored for designing a PVEH connected to a synchronized switch circuit. The governing equation, parameter determination, and numerical analysis technique are useful tools for designing a PVEH connected to a synchronized switch circuit. For simplicity, this study did not consider the voltage drops in the diode rectifier and the dielectric and piezoelectric losses of the piezoelectric materials in modeling. For accurate predictions, these factors should be incorporated into the model depending on the electromechanical coupling states, the magnitude of piezoelectric voltage, and/or acceleration level.

Footnotes

Appendix

Acknowledgements

I thank David MacDonald from Edanz for editing a draft of the manuscript.

Declaration of conflicting interests

The author declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The author received financial support from JSPS KAKENHI [grant number 20K04347].

ORCID iD

Haruhiko Asanuma

References

Asanuma

(2021) Selecting nonlinear piezoelectricity for fully autonomous self-sensing synchronized switch damping on inductor technique. Mechanical Systems and Signal Processing 159: 107846.

Asanuma

Komatsuzaki

(2020) Nonlinear piezoelectricity and damping in partially-covered piezoelectric cantilever with self-sensing synchronized switch damping on inductor circuit. Mechanical Systems and Signal Processing 144: 106867.

Asanuma

Sakamoto

Komatsuzaki

, et al. (2018) Comparative study of electrical and switch-skipping mechanical switch for self-powered SSHI in medium coupled piezoelectric vibration energy harvesters. Smart Materials and Structures 27(7): 075025.

Asanuma

Sakamoto

Komatsuzaki

, et al. (2019) Electromechanical dynamics of Z-shaped mechatronic synchronized switch harvesting on inductor (SSHI) harvester. Japanese Journal of Applied Physics 58(6): 064501.

Asanuma

Yamauchi

(2022) Modal analysis and demonstration of metastructure combining local resonators and an autonomous synchronized switch damping circuit. Journal of Low Frequency Noise Vibration and Active Control 41(2): 540–550.

Badel

Benayad

Lefeuvre

, et al. (2006) Single crystals and nonlinear process for outstanding vibration-powered electrical generators. IEEE Transactions on Ultrasonics Ferroelectrics and Frequency Control 53(4): 673–684.

Brenes

Morel

Gibus

, et al. (2020a) Large-bandwidth piezoelectric energy harvesting with frequency-tuning synchronized electric charge extraction. Sensors and Actuators. A, Physical 302: 111759.

Brenes

Morel

Juillard

, et al. (2020b) Maximum power point of piezoelectric energy harvesters: A review of optimality condition for electrical tuning. Smart Materials and Structures 29(3): 033001.

Chen

Arroyo

, et al. (2017) Shock reliability enhancement for MEMS vibration energy harvesters with nonlinear air damping as a soft stopper. Journal of Micromechanics and Microengineering 27(10): 104003.

10.

Chen

Liu

, et al. (2019) Switching delay in self-powered nonlinear piezoelectric vibration energy harvesting circuit: Mechanisms, effects, and solutions. IEEE Transactions on Power Electronics 34(3): 2427–2440.

11.

Gibus

Gasnier

Morel

, et al. (2022) Non-linear losses study in strongly coupled piezoelectric device for broadband energy harvesting. Mechanical Systems and Signal Processing 165: 108370.

12.

Goldschmidtboeing

Eichhorn

Wischke

, et al. (2011) The influence of ferroelastic hysteresis on mechanically excited PZT cantilever beams. In: 11th international workshop on micro and nanotechnology for power generation and energy conversion applications (PowerMEMS 2011), Seoul, 15–18 November, pp.114–117.

13.

Guyomar

Badel

Lefeuvre

, et al. (2005) Toward energy harvesting using active materials and conversion improvement by nonlinear processing. IEEE Transactions on Ultrasonics Ferroelectrics and Frequency Control 52(4): 584–595.

14.

Hara

Saito

Makihara

(2019) Compact, digital and self-powered piezoelectric vibration energy harvester with generation control using voltage measurement circuit. Sensors and Actuators. A, Physical 299: 111609.

15.

Hara

Yamamoto

Makihara

(2020) Self-sensing state estimation of switch-controlled energy harvesters. Journal of Intelligent Material Systems and Structures 31(20): 2326–2341.

16.

Hegendörfer

Steinmann

Mergheim

(2022) Nonlinear finite element system simulation of piezoelectric vibration-based energy harvesters. Journal of Intelligent Material Systems and Structures 33(10): 1292–1307.

17.

Zhou

Wang

, et al. (2021) Influence of effective electrode coverage on the energy harvesting performance of piezoelectric cantilevers. Energy Conversion and Management 248: 114758.

18.

Ito

Asanuma

Komatsuzaki

(2020) Evaluation of nonlinear parameters in vibration energy harvester based on nonlinear piezoelectricity. In: The Japan society of mechanical engineers (JSME) dynamics, measurement and control division conference on dynamics and design conference, Online, Japan, 1–4 September 2020, paper no. 447. Tokyo: JSME. (in Japanese)

19.

Lallart

Guyomar

(2008) An optimized self-powered switching circuit for non-linear energy harvesting with low voltage output. Smart Materials and Structures 17(3): 035030.

20.

Lallart

Guyomar

(2012) Switching delay effects on nonlinear piezoelectric energy harvesting techniques. IEEE Transactions on Industrial Electronics 59(1): 464–472.

21.

Leadenham

Erturk

(2015) Unified nonlinear electroelastic dynamics of a bimorph piezoelectric cantilever for energy harvesting, sensing, and actuation. Nonlinear Dynamics 79(3): 1727–1743.

22.

Leadenham

Erturk

(2020) Mechanically and electrically nonlinear non-ideal piezoelectric energy harvesting framework with experimental validations. Nonlinear Dynamics 99(1): 625–641.

23.

Lefeuvre

Badel

Brenes

, et al. (2017) Analysis of piezoelectric energy harvesting system with tunable SECE interface. Smart Materials and Structures 26(3): 035065.

24.

Lefeuvre

Badel

Richard

, et al. (2005) Piezoelectric energy harvesting device optimization by synchronous electric charge extraction. Journal of Intelligent Material Systems and Structures 16(10): 865–876.

25.

Wang

Cui

, et al. (2022) Recent progress and development of interface integrated circuits for piezoelectric energy harvesting. Nano Energy 94: 106938.

26.

Liang

Liao

(2012) Improved design and analysis of self-powered synchronized switch interface circuit for piezoelectric energy harvesting systems. IEEE Transactions on Industrial Electronics 59(4): 1950–1960.

27.

Machado

Febbo

Gatti

, et al. (2020) A piezoelectric beam model with geometric, material and damping nonlinearities for energy harvesting. Smart Materials and Structures 29(9): 095009.

28.

Mahmoodi

Jalili

Daqaq

(2008) Modeling, nonlinear dynamics, and identification of a piezoelectrically actuated microcantilever sensor. IEEE/ASME Transactions on Mechatronics 13(1): 58–65.

29.

Mak

McWilliam

Popov

, et al. (2011) Performance of a cantilever piezoelectric energy harvester impacting a bump stop. Journal of Sound and Vibration 330(25): 6184–6202.

30.

Makihara

Asahina

(2017) Analog self-powered harvester achieving switching pause control to increase harvested energy. Smart Materials and Structures 26(5): 055007.

31.

Makihara

Yamamoto

Yoshimizu

, et al. (2015) A novel controller to increase harvested energy from negating vibration-suppression effect. Smart Materials and Structures 24(3): 037005.

32.

Meesala

Hajj

(2018) Identification of nonlinear piezoelectric coefficients. Journal of Applied Physics 124(6): 065112.

33.

Morel

Badel

Wanderoild

, et al. (2018) A unified N-SECE strategy for highly coupled piezoelectric energy scavengers. Smart Materials and Structures 27(8): 084002.

34.

Poblete

Peralta

Ruiz

(2021) Tuning nonlinear model parameters in piezoelectric energy harvesters to match experimental data. ASCE-ASME Journal of Risk and Uncertainty in Engineering Systems, Part B: Mechanical Engineering 7(1): 010904.

35.

Richard

Guyomar

Lefeuvre

(2007) Self-powered electronic breaker with automatic switching by detecting maxima or minima of potential difference between its power electrodes. Patent, publication number: WO/2007/063194.

36.

Shivashankar

Kandagal

(2019) Characterization of elastic and electromechanical nonlinearities in piezoceramic plate actuators from vibrations of a piezoelectric-beam. Mechanical Systems and Signal Processing 116: 624–640.

37.

Shu

Lien

(2007) An improved analysis of the SSHI interface in piezoelectric energy harvesting. Smart Materials and Structures 16(6): 2253–2264.

38.

Stanton

Erturk

Mann

, et al. (2010a) Nonlinear piezoelectricity in electroelastic energy harvesters: Modeling and experimental identification. Journal of Applied Physics 108(7): 074903.

39.

Stanton

Erturk

Mann

, et al. (2010b) Resonant manifestation of intrinsic nonlinearity within electroelastic micropower generators. Applied Physics Letters 97(25): 254101.

40.

Stanton

Erturk

Mann

, et al. (2012) Nonlinear nonconservative behavior and modeling of piezoelectric energy harvesters including proof mass effects. Journal of Intelligent Material Systems and Structures 23(2): 183–199.

41.

Wagner

Hagedorn

(2002) Piezo-beam systems subjected to weak electric field: Experiments and modelling of non-linearities. Journal of Sound and Vibration 256(5): 861–872.

42.

Wang

Gao

, et al. (2020) Nonlinear analysis of piezoelectric wind energy harvesters with different geometrical shapes. Archive of Applied Mechanics 90(4): 721–736.

43.

Weisstein

(2022) Quartic equation. Available at: mathworld.wolfram.com/QuarticEquation.html (accessed 3 June 2022).

44.

Yang

Upadrashta

(2016) Modeling of geometric, material and damping nonlinearities in piezoelectric energy harvesters. Nonlinear Dynamics 84(4): 2487–2504.

45.

Yang

Tang

, et al. (2020) Modelling and analysis of an out-of-plane electret-based vibration energy harvester with AC and DC circuits. Mechanical Systems and Signal Processing 140: 106660.

46.

Zhang

Xiang

Tang

(2021) Modeling, analysis and comparison of four charging interface circuits for piezoelectric energy harvesting. Mechanical Systems and Signal Processing 152: 107476.

47.

Zhao

Yang

(2015) Analytical solutions for galloping-based piezoelectric energy harvesters with various interfacing circuits. Smart Materials and Structures 24(7): 075023.

48.

Zhao

Yang

(2017) Comparison of four electrical interfacing circuits in wind energy harvesting. Sensors and Actuators. A, Physical 261: 117–129.