Utilization of Stochastic P-Bifurcation for Simultaneous Energy Harvesting and Vibration Suppression: An Experimental Investigation

Bi-stable energy harvesters have attracted significant attention over the last one and a half decades [1–3]. This is especially the case for systems with random broadband excitation as a larger bandwidth of operation can improve efficiency. Enhanced energy harvesting due to a jump between potential wells (interwell oscillation) [4–6] as well as performance degradation in linear suspension systems when purposed to harvest energy [7–9] have marked a couple reasons for further exploration of bi-stable systems.

The experimental exploration of bi-stable energy harvesters subjected to random excitation has accompanied much of the research in this particular field. For example, in Ref. [6], the authors qualitatively verified simulation results for a Duffing-type harvester subjected to band-limited noise excitation, investigating how the center frequency and the bandwidth of the excitation affected the output voltage in the mono-stable and bi-stable configurations. In Ref. [4], the researchers experimentally validated simulation results for a bi-stable energy harvester under band-limited white noise excitation which was transformed into a flexible bi-stable energy harvester exhibiting a variable potential energy function due to a magnet on a flexible beam. Similarly, in Ref. [3], the authors experimentally validated simulation results for a bi-stable two degrees-of-freedom harvester system and investigated the effect of varying excitation intensity, bandwidth, and center frequency. Likewise, in Ref. [5], the authors explored and experimentally validated the theoretical results for a piezoelectric bi-stable harvester subjected to wide-band random excitation, while varying excitation intensity and the separation between equilibrium positions. In short, Refs. [3–6] and others like Ref. [10] all share a common conclusion: bi-stable energy harvesting systems will outperform mono-stable counterparts under broadband random excitation precisely when the excitation is large enough to produce interwell oscillations and they all experimentally validate or verify their theoretical results. Nevertheless, prediction of the parameter values required for interwell oscillations has remained open thus far.

When interwell oscillations are permitted, the multi-stability in energy harvesters in the case of random excitation manifests itself as a multimodal (multiple peaks) probability density function (PDF) associated with the states of the system [3,11,12]. As such, analysis often involves obtaining the PDF by solving the Fokker–Planck–Kolmogorov equation (FPE) with the finite element method [12,13] or by obtaining an approximation with the Wiener path integration method (WPI) [14–17]. Otherwise, statistical moment equations governing the system’s statistics can also be derived [18,19] as another means of analysis. A Monte Carlo simulation (MCS) is then often used for verification purposes [17,19–21], which allows for the determination of the necessary statistics with direct numerical simulation.

To understand when multiple peaks in the PDF can occur, or equivalently when interwell motion capabilities change for multi-stable systems, one can understand P-bifurcations. Literature associated with stochastic dynamics such as Ref. [22], which provides an overview of the advancements in stochastic dynamics by the close of the twentieth century, defines P-bifurcation as qualitative changes in the stationary PDF and the related Markov process. For a more in-depth understanding of the mathematical methods used to examine P-bifurcations, one can refer to works such as Ref. [23], which compiles the stochastic dynamics literature from the late twentieth century. A notable article in the text is one on P-bifurcation in the noisy Duffing–van der Pol equation [24]. The authors of this work noted that the random motion was fast in nature and occurred as random oscillations or rotations along the trajectories of the deterministic system. They further used this fact to justify the use of stochastic averaging [25] and then proceeded to derive the stationary PDF from the FPE. This allowed them to mathematically determine when the extrema in the PDF changed, indicating P-bifurcation. On the other hand, the authors in Ref. [26] begin by reviewing deterministic center manifold theory, utilizing unfolding parameters to examine behavior around bifurcations such as Hopf, transcritical, and pitchfork types. They then derive the FPE and its solution near the center manifold to identify P-bifurcations. Analyses cover both additive noise, where noise is added to the state space, and multiplicative noise, where noise multiplies the states in the system. In addition to the aforementioned analysis methods, the authors of Ref. [27] mentioned the notion of bifurcations in the moment equations governing the statistical moments associated with the system [25,28]. The method employed in the current work, similar to that in Ref. [19], involves detecting the change in stability of moment equation fixed points by determining when the determinant of the Jacobian vanishes.

Some approaches in the past to study bifurcations have included the numerical computation of the PDF with the finite element method [13] or deriving estimates for the PDF [13,29–31], which have been computationally limited to systems with less than five states. Then in Ref. [32], the authors used a method based on Shannon Entropy while Cosner and Tai used a WPI-based PDF approximation algorithm for a system of seven states in Ref. [17]. To the authors’ knowledge, no openly available research attempts to predict P-bifurcation pertaining to experiments with wide-band random excitation. Lastly, while many researchers have remedied interwell oscillation difficulties by studying tri-stable [33,34] or even quad-stable [35] energy harvesters to boost energy harvesting efficiency, it is hypothesized from evidence in Refs. [17,19,21] that by using inertial nonlinearity rather than stiffness nonlinearity, a bimodal PDF could still be produced without a multi-stable potential and is often associated with large amplitudes of oscillation just as in multi-stable systems with interwell motion. The current work will exploit such a system and seek to predict changes in the number of peaks in the PDF in experiments, indicating when the mentioned large amplitude oscillations may occur.

Recently, Cosner and Tai studied a device known as an inerter pendulum vibration absorber (IPVA), which can exhibit bi-modality with inertial nonlinearity. Through a theoretical and numerical investigation, it was found that this device could effectively reduce vibrations in a single-degree-of-freedom structure subjected to Gaussian white noise base excitation [19]. It was determined that P-bifurcation associated with the marginal PDF for the pendulum angle of oscillation was directly related to optimal vibration suppression. In Ref. [19], they also developed an algorithm which tracked P-bifurcation of the device through determining the singularity of the Jacobian of the statistical moment equations [25,28]. The same device was connected to a generator and served as an energy harvesting shock absorber in a two degrees-of-freedom quarter car model in Ref. [21]. In this case, it was used for simultaneous optimization of energy harvesting and passenger ride comfort which amounted to a 45$%$ improvement over the linear benchmark system before adding any model predictive control. The work in Refs. [19,21] served as motivation for further exploration of the device by Cosner and Tai in Ref. [17]. In Ref. [17], they predicted P-bifurcation of the quarter car model in Ref. [21] using a WPI algorithm uniquely adapted from work in Ref. [14], and correlated the P-bifurcation with great performance in terms of road handling, ride comfort, and energy harvesting objectives. In this scenario, the implementation and tuning of the IPVA near the P-bifurcation point led to enhancements of 40% in power, 60% in ride comfort, and 60% in road handling. The study further clarified the power improvements by analyzing the power spectral density (PSD) of the generator’s velocity, which showed a bandwidth increase and thus enhanced power output when the system states’ PDF exhibited a bimodal distribution. In Ref. [21], the increased bandwidth in the PSD was shown to be coupled with the annihilation or reduction of peaks in the PSD associated with sprung mass acceleration, resulting in better ride comfort.

The research surrounding the multi-objective performance capabilities of the nonlinear IPVA when subjected to random excitation such as in Refs. [17,21] is very relevant to modern research. As pointed out in a recent review of dual-objective harvesters [36], the research surrounding the use of nonlinearity to achieve multi-objective performance is still scarce. It is even more scarce when considering excitation which is stochastic in nature. The research that does exist is often associated with the longtime running research topic of nonlinear energy sinks (NES) [37] which inherently require cubic stiffness nonlinearity via the use of elastic components or magnetic arrangements. For example, the researchers in Ref. [38] theoretically coupled an NES with a giant magnetostrictive energy harvester and applied the device to a single-degree-of-freedom structure subjected to Gaussian white noise excitation to achieve energy harvesting and vibration suppression objectives simultaneously. Regardless, the authors of the current work wish to avoid stiffness nonlinearities which could make the system too bulky and impractical if applied to an automobile suspension system.

The focus of this work will be on the bimodal IPVA energy harvesting device with an application to a single-degree-of-freedom structure subjected to Gaussian broadband excitation. The objective is to conduct an experimental investigation into P-bifurcation of the device and the correlation with simultaneous energy harvesting and vibration suppression. Notably, this study will include the experimental parameter characterization for a developed prototype to accommodate accurate predictions of energy harvesting potential, vibration suppression, and bifurcation boundaries. Adequate experimental data will be used to verify predictive capabilities, confirm the bifurcation phenomenon, and gather information regarding power harvested and vibration suppression.

The rest of this work is organized as follows. In Sec. 2, the IPVA device, system model, and equations of motion (EOM) are fully introduced. In Sec. 3, the experimental setup is introduced and a characterization is done to determine unknown parameters and allow for accuracy in predictions. In Sec. 4, the moment equations are introduced and followed by an analytically determined bifurcation diagram in the noise intensity—electrical damping plane with the help of the tracking algorithm developed in Ref. [19]. Additionally, experimentally determined system parameters are used with the bifurcation diagram to predict bifurcation in the experiment. The harvested electric power and vibration mitigation are quantified and compared to a linear benchmark system without pendulums, while correlations are also made between the PDFs and performance. The experimental PDFs are ultimately compared with an MCS. The work is concluded in Sec. 5.

In this section, the IPVA device and corresponding system model are introduced. The model for the system is shown in Fig. 1. The system has a suspended mass $Ms$ with a degrees-of-freedom $xs$⁠, stiffness $ks$⁠, mechanical damping $cm$⁠, and an attached IPVA device with planetary gear coupling to a DC generator. The system also has a base motion $xr$⁠.

A holonomic constraint is imposed by the ball screw giving $xs−xr=Rθ$⁠, where R is the lead value of the screw divided by $2π$⁠. Figures 1(b) and 1(c) then show the bottom and top of a carrier that holds two pendulums, each of mass $mp$ and principal inertia $Jp$ about their centers of mass. Additionally, the center of mass of each pendulum is displaced at a distance $r$ from the respective pivot points which are at a distance $Rp$ from the center of the carrier. These pendulums are free to rotate with an angle $ψ$ relative to the carrier, while the carrier with principal moment of inertia $Jc$ about its vertical axis is free to rotate with an angular degrees-of-freedom $θ$ about the ball screw’s vertical axis. Also attached to each pendulum is then a planetary gear of mass $mgear$ and principal moment of inertia $Jpg$⁠. These planetary gears are then meshed with a sun gear connected to the input of a gearbox whose output is connected to the input of a DC generator. The gearbox shaft has a degrees-of-freedom $ψg=θ−ψ$ due to the holonomic constraint imposed by the planetary gear system. The motor rotor then has a degrees-of-freedom $ngψg$⁠, where $ng$ is the gear ratio. The principal moment of inertia of the sun gear and generator rotor are $Jsg$ and $Jr$⁠, respectively.

In Eqs. (3) and (4), $x=(θ,ψ)T$⁠, $M(ψ)$ is the dimensionless inertia matrix with nonlinear dependence on $ψ$⁠, $C$ is the dimensionless damping matrix, $K$ is the dimensionless stiffness matrix, and the vector $f$ is associated with stochastic excitation. The excitation term $W(τ)=−x¨rMs/(MsRω02)$ is then defined to be the normalized Gaussian white noise with zero mean and constant two-sided power spectral density $D$⁠. Here $D=d/(ω03R2)$⁠, where $d$ is the physical noise intensity in units of squared acceleration per frequency. Furthermore, $g$ represents the nonlinear Coriolis and centrifugal terms. It is finally apparent that $μr$ and $η$ are the quantifiers of nonlinearity in the system.

The experimental setup correlating with Fig. 1 is now shown in Fig. 2. The setup consists of the IPVA, a generator, a mass suspended with eight coil springs and clamps, a signal conditioner, an electrodynamic shaker, an amplifier, a vibration controller, and a pendulum tracking system. In this section, the measured or experimentally determined values for all parameters outlined in the previous section are given.

With the experimental setup shown in Fig. 2, a three-hour random excitation experiment was done for the case of four different load resistances attached to the generator, giving four different electric damping (⁠$ξe$⁠) values with the use of Eqs. (1) and (2). This is outlined in Table 1. The resistor values were chosen to incorporate a sufficiently large electrical damping range while using resistors that were immediately available. Note that open circuit implies infinite load resistance. Also, note that a total of 30 min of data was removed in order to eliminate transient beginning and ending dynamics which might have a misleading effect when interpreting some statistical information. The total duration of the analyzed data in each case was then about 2.5 h.

The excitation used in the experiment was such that the base acceleration had a root mean square (RMS) value of $0.8g$ where g is the acceleration due to gravity and the corresponding PSD is shown in Fig. 3. Note that this PSD was constructed to be broadband on the frequency band between 2 and 13 cycles per second while also including the resonant frequency determined in Ref. [39] and also within the voltage, current, and stroke limits of the electrodynamic shaker. Shaker current, voltage, and minimum frequency limitations are directly proportional to the acceleration, velocity, and displacement output limitations, respectively.

For the purpose of verifying bifurcation predictions, the same experiments were also done for the case of $0.6g$ RMS base acceleration as well as a $0.8g$ RMS base acceleration case with a shorted circuit (zero load resistance). The corresponding PDFs will be analyzed and it shall be determined whether the bifurcation boundaries can be used to predict their qualitative nature.

The linear system frequency response function (FRF), corresponding to the system in Fig. 2 with pendulums removed, was accurately fitted with experimental data to obtain the values of $ω0$⁠, $ξ$⁠, and $α$⁠. This utilized a sine sweep experiment and is elaborated in the previous work by Ref. [39]. Various fitted, measured, and generator parameters used in this study are presented in Tables 2 and 3. It should be noted that $μr$ was directly measured, while $μp$ and $η$ were computed using solidworks. Additionally, $ce$ varies with load resistance as described by Eq. (2), but once a load resistance is selected, it is calculated using the reported values for the Maxon 110207 4.5 W DC motor with a gear ratio of $ng=5.2$⁠. Finally, although mechanical damping $ξ$ was determined through the FRF fit, it was observed that it significantly differs from the value obtained in the actual experiment described in Sec. 3.1; hence, it is discarded hereinafter. It should finally be noted that stiffness $ks$ and mass $Ms$ in $ω0$ were chosen to keep the linear resonant frequency low and accommodate a lower excitation frequency upper bound while also keeping the mechanical damping ratio $ξm$ as low as possible with the springs and mass which were available on hand. Pendulum parameters $μr$ and $η$ were chosen after numerical simulation predicted that they were sufficiently large to accommodate the bimodal phenomenon of interest while also remaining relatively compact. The choice of generator torque constant $κt$⁠, internal resistance $Ri$⁠, and gear ratio $ng$ was governed by a combination of what was commercially available and the intent to limit generator inertia $μg$⁠, limit internal resistance, and use a gear ratio which would put the generator speed near its nominal rated value in the documentation.

To execute the optimization described in Eq. (5), an initial three-dimensional grid of finite parameter values was established. Initial fitting in the time domain detailed in Ref. [39] helped narrow down the parameter space to the range specified in Eq. (5). Next, parallel computing leveraging high-performance computing resources with $4×128$ CPU cores and 921.6 GB of memory per 128 cores to solve the EOM in Eq. (5) with ode45 in matlab. Consequentially, the objective function value was determined. The minimal objective function that satisfied the constraint was then determined and the optimized damping parameters for all resistant loads along with the corresponding maximum RMS velocity percent errors are presented in Table 4. It should be noted that each resistance load case was associated with a unique optimal $Ξ$⁠. Throughout the rest of this study, $Ξ$ is considered as a function of $ξe$ (dependent on load resistance). This approach is necessary because the model described in Eq. (5) does not fully capture the complex dynamics involved. For instance, intricate dynamics associated with ball screw damping, discussed in Refs. [41,42], are omitted.

Following the optimization routine of Eq. (5), a comparison between experimental and simulated PSDs is displayed in Fig. 4. Note that the maximum RMS velocity error is below 3% across all PSDs and the low frequency contributions are accurate qualitatively, while the peak PSD values do not align as closely. It was proven difficult to obtain complete agreement. While this is deemed acceptable given the primary focus on the RMS error in this research, it highlights some flaws in the mathematical model.

In concluding this section, it is noted that all RMS velocities and mean power, with power being directly related to the square of the generator velocity [17,21], were within a 3% error margin, which the authors consider to be sufficiently accurate. Consequently, it is anticipated that predictions regarding vibration suppression (reduction in relative velocity of $Ms$⁠) and energy harvesting will be significantly accurate, at least when electrical damping and noise intensity are relatively close to the experimental values in magnitude. Additionally, it is important to highlight that the optimization that initially only minimized RMS velocity error resulted in RMS error as low as 0.62%. However, predictions related to the rotational dynamics were poor. Hence, the authors have opted to use the objective function in Eq. (5) for enhanced accuracy in predicting complicated dynamics and improving the likelihood of forecasting P-bifurcation in the following section.

Cosner and Tai [19] have already formulated a tracking algorithm that can be used to predict the critical values in an arbitrary two-parameter plane for which a P-bifurcation, or change between the total number of PDF maxima, may occur. In this section, the stochastic bifurcation analysis is briefly introduced and the prediction is compared with experimental results.

To alleviate the infinite hierarchy of dependence of moments on higher order moments, as detailed in Ref. [19], $n$ is set equal to two, and Gaussian closure is used [25,28]. Next, the Jacobian of the moment equations in Eq. (6) with respect to moment variables and $ϕ0$⁠, is evaluated at $ϕ0=0$ and the parameter values required for a vanishing Jacobian determinant are determined. Then with the algorithm in Ref. [19] all P-bifurcation solution candidates are readily found.

The dimensionless value is then $D=0.2687$⁠.

Using the noise intensity, $ξe$⁠, $Ξ(ξe)$ from Table 4, and the bifurcation tracking algorithm in Ref. [19], a bifurcation diagram was generated for each resistor case in Table 1 with the RMS acceleration excitation of $0.8g$⁠. The maximum and minimum boundary for a given electrical damping was then determined and the resultant boundary region is shown in Fig. 5. In this figure are also markers indicating the experimental damping and noise intensity values, including the values corresponding to the experiments only used to test the boundary accuracy. It is assumed that the region between the curves in Fig. 5 serves as an estimate for where p-bifurcation can occur. In other words, an electrical damping-noise intensity pair outside of this region is assumed to be not correlated with a bifurcation point.

The marginal PDFs for the pendulum angle corresponding to each experimental parameter set were also obtained via a Monte Carlo simulation with $∼4⋅104$ realizations, a dimensionless simulation time of $104$⁠, and a time-step of $10−3$⁠. Note that this simulation time corresponds to only $∼130s$ in physical time which is drastically smaller than the experimental time of 2.5 h. As a result, any dynamics that develop on a timescale larger than 130 s might not be accurately captured. However, given the number of realizations needed to reach stationarity, running the MCS with a simulation time of 2.5 h is impractical within a reasonable timeframe. With this disclaimer aside, the resultant marginal PDFs are shown in Fig. 6 and the corresponding marginal PDFs from the experiment are shown in Fig. 7. Note that the pendulum will undergo oscillations centered around $ψ=2πk$ rad for any integer k with possible intermittent rotations in between. As such, the pendulum angle $ψ$ has been replaced in the PDFs with $ψmod2π$⁠. This does have the side effect of the separation of the PDFs from the $ψ$ axis near $ψ=±π$ when rotation exists.

Initially, it is important to recognize that the bifurcation diagram in Fig. 5 delineates both right and left boundary regions. Considering the relationship between damping and resistance in Table 1, transitioning from a $3.9Ω$ to a $1Ω$ load resistor necessitates entering the right bifurcation boundary region from the left and ending near the outer edge. With that being said, the bifurcation boundary generated with the parameters from the $1Ω$ load resistor case actually had an outer edge that was to the left of the experimental data point. In this case, one can assume that bifurcation has in fact occurred. Notably, the PDF in Fig. 6 associated with the $1Ω$ case has negative curvature at its origin and is distinctly monomodal, while the PDF for the $3.9Ω$ case on the left side of the boundary region shows positive curvature at the origin and is distinctly bimodal. The experimental PDFs in Fig. 7 exhibit similar characteristics in terms of curvature at the origin, with more pronounced bi-modality. Furthermore, the experimental short circuit case is shown in Fig. 7, corresponding to the marker with the highest damping value in Fig. 5. The PDF is indeed more distinctly monomodal, indicating that the transition has fully occurred from bimodal to monomodal. The right bifurcation boundary has been predicted reasonably well. Moreover, if the left boundary region is accurate, the $10Ω$ case is likely to display a bimodal PDF as well. Both simulation and experimental PDFs in Figs. 6 and 7 reveal this bi-modality, also with rotation indicated by the separation from the $ψ$-axis at $ψ=±π$⁠, attributed to reduced electrical damping. Crossing the left boundary region to examine the open circuit scenario, one might anticipate another qualitative shift in the PDFs. In this instance, Figs. 6 and 7 suggest a possible flattening of the PDFs, while they show an apparent change from positive to negative curvature at the origin, indicating bifurcation. The precise nature of this change remains uncertain, yet a comparable qualitative change was noted in the authors’ recent publication [17], where the left boundary was presumed to lead to pure rotation. Lastly, Fig. 8 shows the experimental PDFs for $0.6g$ RMS acceleration. The boundary would predict that all of these should be bimodal, with the exception of the open circuit case which lies in the region of uncertainty. After closer examination, the three bimodal predictions are in fact verified, albeit more asymmetric due to frictional effects or not enough time to become stationary. The PDF with a prediction in the bifurcation region of uncertainty then appears be possibly locally bimodal with a somewhat flattening of the PDF. This indicates a potential transition toward pure rotation.

Some concluding observations pertaining to this subsection are now provided. The PDFs derived from the MCS were produced using numerous realizations in the time domain, in contrast to the PDFs associated with the experiment, which were derived from a single 2.5 h experiment. Consequently, it is highly probable that the PDFs depicted in Fig. 6 are closer to being stationary compared to those in Fig. 7, as suggested by the greater symmetry in Fig. 6. However, there remains a possibility that the duration of the simulation was insufficient to achieve fully developed PDFs, even if they are relatively symmetric compared to the experimental PDFs. Furthermore, it is important to highlight that generating these PDFs took around 4 h, leveraging high-performance computing resources with $8×128$ CPU cores and 921.6 GB of memory per 128 cores. Additionally, it should be noted that the Gaussian white noise excitation utilized in the analysis is merely an approximation of the band-limited noise employed in the experiment. Lastly, the precision of the pendulum tracking system was constrained by a sampling rate of five hundred samples per second, yet its measurement accuracy was confirmed to be within a couple of degrees. Nevertheless, the qualitative change in the PDF using the approximate bifurcation boundary was predicted reasonably well.

In this section, the PSDs are presented for the relative velocity of the suspended mass and harvested power computed with the data obtained from the experiments discussed in Sec. 3.2.1. Additionally, the same experiments were done for a linear benchmark system where the pendulums are removed, but equivalent inertia is added so that the natural frequency remains about the same for a fair comparison. However, note that the pendulums act as vibration absorbers and so their removal results in very high displacement amplitudes of the suspended mass. For safety and to preserve the integrity of the experiment, the linear system was excited at only half of the root mean square acceleration excitation and the resultant amplitude was scaled by two. This is simply utilizing the principles of linearity, i.e., a force twice as large will give a response twice as large. Ultimately, after this scaling, the PSDs from the linear system experiment will be compared to those from the experiment with the addition of the IPVA.

First, Fig. 9 shows the PSD for the electric power harvested across the resistant load. Note that the power in this work was considered as the power across the resistive load which was connected in between the generator terminals. The power harvested was calculated directly from measured voltage $(Vm)$ as $P=Vm2/Rl$⁠. In this scenario, it becomes evident that the $1Ω$ resistance setup correlates with the least electric power near the resonant frequency, while it still outperforms all linear system cases in terms of mean power, calculated as the area under the PSD. A closer look at the PDF in Fig. 7 reveals that this $1Ω$ resistance setup is also a monomodal configuration. The correlation between a monomodal PDF and suboptimal energy harvesting exists. For $Rl=10Ω$⁠, the mean power is around four times that of the best linear system, and the power near the resonant frequency is more than double that of the best linear system (⁠$Rl=3.9Ω$⁠). A significant portion of the mean power can be attributed to the contributions from low frequency components in the PSD below the excitation frequency range due to the pendulums’ rotation and other nonlinear effects. Additionally, it should be noted that the actual amount of power is not discussed in this work as the concern is primarily with how the power correlates with the characteristics of the PDF and how the nonlinear system performance compares to that of the linear benchmark system.

Using the parameters from Table 3, one finds when rounded to the hundredths place, $Rlopt=3.17Ω$⁠, which is very close to the experimental value of $3.9Ω$⁠.

It is clear in Fig. 10 that the nonlinear system once again outperforms in terms of vibration suppression near the resonant frequency except when there is an open circuit condition and electrical damping is zero. A reduction by around a factor of four in the relative velocity PSD value is seen near the resonant frequency and an overall halving in the mean squared value of the relative velocity has been computed for the case of $Rl=1Ω$ in the system with the IPVA. As long as power harvesting is a priority and the circuit is not open, it appears as though the IPVA will allow for superior vibration suppression and power harvesting over the linear system. This possibility for simultaneous objectives somewhat echoes the result of utilizing the device in a two degrees-of-freedom quarter car model such as in Refs. [17,21].

In concluding this section, it is important to solidify the hypothesis for why energy harvesting appears to be optimal when the nonlinear system PDF is bimodal. More specifically, the most superior power was found in the $10Ω$ case. A closer look at Figs. 6 and 7 reveals that this case was also associated with a smaller decrease in PDF magnitude in between the peaks when compared to the bimodal PDF associated with the $3.9Ω$ resistance case. It is hypothesized that this smaller decrease in PDF magnitude between peaks is responsible for more frequent large amplitude jumps between the pendulum angles specified by the PDF extrema. The work associated with the same energy harvesting IPVA device implemented in a quarter car suspension system also concluded that energy harvesting performance was generally better when near bifurcation, where the PDF magnitude reduction is minimal [17]. Furthermore, this hypothesis would be analogous to how a shallow well in double potential-based systems allows for enhanced interwell oscillations which is known and exploited by many researchers such as in Ref. [10].

In this study, an inertially nonlinear device was equipped with a DC generator to conduct experiments on and verify predictions regarding the modality of the marginal PDF of the pendulum angle of oscillation as well as the device’s ability to suppress vibrations and harvest energy when applied to a single-degree-of-freedom structure under Gaussian broadband base excitation, treated analytically as white noise. The experimental excitation, closely approximating Gaussian white noise within experimental constraints, involved a broadband excitation with a root mean square acceleration of $0.8g$ and a bandwidth of eleven cycles per second, incorporating the system’s resonant frequency and a minimum value of two cycles per second.

To facilitate experimental verification, the system’s unknown parameters were formally characterized with frequency domain optimization, where the optimization routine was designed to minimize the mean squared error in the pendulum velocity on the frequency band below two cycles per second while constraining the RMS velocity discrepancy between the simulations and actual experiments to be below 3%. All parameters whether measured, calculated, or characterized via the optimization routine were documented. Experiments utilized four different load resistance settings: open circuit, $1Ω$⁠, $3.9Ω$⁠, and $10Ω$⁠. Each resistance setting was associated with a unique set of damping parameters determined from optimization. The RMS velocities and the corresponding mean power (proportional to the squared RMS velocity) were shown to be predicted within a margin of error of 3%. This signifies accurate estimations of both power generation and vibration reduction. Additionally, a visual comparison between the simulation and the PSDs for all velocities demonstrated excellent agreement.

A numerically determined P-bifurcation boundary in the plane of noise intensity versus electrical damping was derived for each experiment at $0.8g$ RMS acceleration excitation, excluding an additional test scenario with a short circuit, and a bifurcation region was established based on all four boundaries. Two clear qualitative changes were predicted and experimentally verified when noise intensity remained constant and electrical damping varied. The right boundary region indicated a transition from a bimodal to monomodal PDF, while the left boundary region suggested a potential flattening of the PDF. The validity of the bifurcation boundary was enhanced with the addition of the short-circuited experiment at $0.8g$ and four additional $0.6g$ experiments. There was significant qualitative agreement between experimental and simulation PDFs.

The power spectral density for electric power in the IPVA system was compared to a linear benchmark system. Results showed that power harvested by the IPVA increased by a factor of two near the resonant frequency, and mean power by a factor of four, compared to the best linear system tested. Additionally, an analytically optimized linear system was shown to still produce less than one third of the power harvested with the IPVA system. The IPVA configuration with the lowest performance in terms of power harvesting had the lowest load resistance of $1Ω$⁠, associated with a monomodal PDF.

Finally, the power spectral density for the relative velocity of the suspended mass was assessed and compared to the most effective experimentally tested linear system for power harvesting, which had a load resistance of $3.9Ω$⁠. All configurations except the open circuit configuration outperformed the linear system, by as much as a factor of four at the peak PSD value. Then all except the open circuit and $10Ω$ resistance configuration had reduced mean square relative velocity, by as much as a factor of two. It was determined that the $10Ω$ resistance configuration had 6% higher mean square velocity than the tested linear system with $3.9Ω$ resistive load. However, further analysis showed that a linear system that was fully optimized for power had 20% higher mean square velocity and still less than one third of the power when compared to the $10Ω$ resistance IPVA configuration.

This material is based upon work supported by the Office of Naval Research under Award No. N00014-22-1-2533. Any opinions, findings, conclusions, or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the funding agency.

• Office of Naval Research (Award No. N00014-22-1-2533; Funder ID: 10.13039/100000006).

There are no conflicts of interest.

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.