Uncertainty Quantification of Piezoelectric Energy Harvesters from Aeroelastic Vibrations

A stochastic approach is presented to evaluate the uncertainties associated with variations in design parameters of a piezoaeroelastic energy harvester. The sensitivities of the harvested power to variations in the load resistance, the eccentricity (distance between the center of mass and the elastic axis), and the nonlinear coefficients are also determined. Moreover, the non-intrusive formulation of the polynomial chaos expansion in terms of the multivariate Hermite polynomials was employed to quantify the sensitivities in the harvested power and the plunge and pitch motions. The results show that the relationship between the input parameters and the harvested power is highly nonlinear. The results show also that the generated power is most sensitive to variations in the eccentricity and that the nonlinear coefficient of the plunge spring is less influential than the nonlinear coefficient of the torsional spring on the harvester’s performance.


Introduction
The use of aeroelastic vibrations to harvest energy has been the topic of recent studies.Bryant and Garcia [1] demonstrated experimentally that energy can be harvested from aeroelastic vibrations using an airfoil section attached to a cantilever.Erturk et al. [2] used a lumped-parameter wingsection model to determine the effect of piezoelectric power generation on the linear flutter speed.Abdelkefi et al. [3] considered the nonlinear behavior of piezoaeroelastic energy harvesters and focused on avoiding subcritical (unsafe) Hopf bifurcations.
The level of the harvested power is a function of the harvester parameters and wind speed.If the energy is harvested from the plunge motion, it would be important to design an aeroelastic system where the amplitude of the pitch motion is kept small and most the wind energy is picked up by the plunge motion.Because the amplitudes of these motions depend on many system parameters, it is important to evaluate the uncertainties associated with variations in the different system parameters and the respective sensitivities to these variations.Although one could independently vary one parameter and keep the others constant to achieve this objective, a better analysis would allow for combined variations in the system parameters.To this end, we quantify the uncertainties associated with the performance of energy harvesters and the sensitivity of the generated power to variations in the harvester's parameters.Particularly, the sensitivities of the harvested power to the load resistance, the eccentricity (distance between the center of mass and the elastic axis), and the nonlinear coefficients of the springs are evaluated.a e-mail: mhajj@vt.eduFig. 1.Schematic of a piezoaeroelastic system under uniform airflow.

Representation and modeling of the energy harvester
The energy harvester, considered in this work, consists of a rigid airfoil that is allowed to move with two degrees of freedom, as shown in Figure 1.The airfoil is supported by linear and nonlinear torsional and flexural springs with a piezoelectric coupling attached to the plunge degree of freedom.The governing equations of this system are written as [2,3]: where m T is the total mass of the wing, including its support structure; m W is the wing mass alone; I α is the mass moment of inertia about the elastic axis; b is the half chord length; x α is the dimensionless distance between the center of mass and the elastic axis; c h and c α are, respectively, the plunge and pitch structural damping coefficients; L and M are the aerodynamic lift and moment about the elastic axis; R is the load resistance; V is the voltage across this and k α0 are the linear structural stiffnesses for the plunge and pitch degrees of freedom, respectively; and k h2 and k α2 are, respectively, the nonlinear structural stiffnesses for the plunge and pitch degrees of freedom.The aerodynamic loads are evaluated using a quasi-steady approximation with a stall model and written as where U is the freestream velocity, c lα and c mα are the aerodynamic lift and moment coefficients, and c s is a nonlinear parameter associated with stall.The effective angle of attack due to the instantaneous motion of the airfoil is given by where a is a constant.To determine the stability characteristics, we introduce the following state variables: The equations of motion are then rewritten in the following vector form [3]: where C(X, X, X) is a cubic vector function of the state variables and The matrix A(U) has a set of five eigenvalues λ i , i = 1, 2..., 5, which determine the stability of the trivial solution of equation (7).We arrange these eigenvalues so that the first four are complex conjugates (λ 2 = λ 1 and λ 4 = λ 3 ).The wind speed at which one or more eigenvalues have zero real parts corresponds to the onset of linear instability and is termed the flutter speed, U f .

Effects of eccentricity and load resistance on the flutter speed
As a first step in the sensitivity analysis, we determine the effects of the eccentricity and the load resistance on the flutter speed.The rest of the parameters of the energy harvester are presented in [4].As shown in Figure 2, the choice of the eccentricity plays an important role in determining the flutter speed.We note that, as the distance between the elastic and gravity axes is decreased to very small values, the flutter speed approaches infinity.The same figure shows also that there is a value for which the flutter speed is minimum.Inspecting Figure 3, we note that the flutter speed depends also on the electrical load resistance.These results are used in the subsequent analysis to make sure that flutter will always take place.To this end, we set the wind speed to 5m/s.

Uncertainties and Sensitivity Analysis of the System Outputs
An important aspect of spectral representation of uncertainty [5] is that it decomposes a random function (or variable) into separate deterministic and stochastic components.
In this representation, any stochastic variable (i.e., α * ), such as the pitch, plunge, voltage, or power is expressed as where α i (t) is deterministic and represents the amplitude of the i th fluctuation, ψ i (ξ) is a random basis function corresponding to the i th mode, and P is the number of output modes.This number is expressed as P + 1 = (n+p)!n!+p! , where p is the order of the polynomial chaos and n is the number of random variables.The choice of the polynomial basis is highly dependent on the distribution of the random variables [5].For Gaussian distribution of the random variables, Hermite polynomials are used as basis functions.A convenient form of the Hermite polynomials is given by where n i=0 c i = p.This set of polynomials form a complete orthogonal set of functions in the random space with respect to the following inner product: (10) with the density function of the n variate standard Gaussian distribution being the weighting function.We note that To determine the polynomial representation of the system output, we need to compute the deterministic coefficients α i .Using the orthogonality condition of the basis functions, we obtain the following relation: Three different methods can be used to obtain approximations of the polynomial expansion coefficients: the sampling based method, the quadrature method, and the linear regression method.The latter allows for more flexibility in the choice of the sampling points and the use of sampling strategies, such as the Latin Hypercube Sampling (LHS).Another advantage of using this method is that the required number of sampling points is very small compared to the other approaches.This is an important advantage particularly for time consuming computations.This method is based on the following matrix equation: which represents the discretized form of equation ( 11).When N is equal to the number P + 1 of polynomials in the expansion, the matrix is square.Therefore, it can be inverted to determine the expansion coefficients α i from the outputs α * i .As aforementioned, we consider uncertainties in R, x α , k α 2 , and k h 2 and study their impacts on the harvested power and pitch and plunge amplitudes.To generate samples, we assume that the uncertainties follow Gaussian distributions.The mean and standard deviation of each considered parameter are given in Table 1.The mean values of the load resistance and the eccentricity are chosen to maintain flutter speeds that are smaller than the considered wind speed in the simulations, which is set equal to 5m/s.Furthermore, the mean values of the nonlinear coefficients of the springs are chosen based on the normal form of the Hopf bifurcation to guarantee a supercritical instability.The α * i are determined numerically by using the governing equations (1-3) for the samples generated based on a Latin Hypercube Sampling methodology.Then, we use the linear regression method to determine the PCE coefficients.Because the relationship between the input and Table 1.The mean and standard deviation values for the load resistance, the eccentricity, and the nonlinear coefficients of the springs ( output parameters is nonlinear, the sensitivity of the output α * to changes in the i-th input parameter ξ i is given by Determining these sensitivities facilitates the identification of the variations that would impact the level of the harvested power.

Results and discussions
Figures 4, 5 and 6 show the empirical density functions of the harvested power, plunge amplitude, and pitch amplitude, respectively and that are associated with the variations presented in Table 2.These probability density functions of the output parameters were computed from 450 realizations based on the Monte Carlo simulations (MCs).We note that, for all system outputs, the seventh-order nonintrusive PCE is needed, as a minimum, to obtain excellent agreement with the MC.This is due to the fact that the relationship between the input parameters and the system outputs is nonlinear.Table 2 shows the sensitivity of the harvested power and the plunge and pitch amplitudes to global variations in the load resistance, the eccentricity, and the nonlinear coefficients of the springs.The results show that increasing the eccentricity and the nonlinear coefficient of the torsional spring results in a decrease in the harvested power and the plunge amplitude.On the other hand, increasing the  load resistance and the nonlinear coefficient of the plunge spring results in an increase in the harvested power.Based on the sensitivity magnitudes in Table 2, global variations in the eccentricity are the most influential on the system outputs.Furthermore, global variations in the load resistance are small in comparison to variations in x α , k α2 , and k h2 for the plunge and pitch amplitudes.In addition, for the chosen mean values, the load resistance has a negligible impact on the harvested power.We note also that the sensitivities of the different outputs to the nonlinear coefficient of the plunge spring are smaller than those of the nonlinear coefficient of the pitch spring and the eccentricity.

Conclusions
We perform a sensitivity analysis of a piezoaeroelastic energy harvester consisting of a pitching and plunging rigid airfoil supported by flexural and torsional springs with a piezoelectric coupling attached to the plunge degree of freedom.We employ the non-intrusive formulation of the polynomial chaos expansion (PCE) in terms of the multivariate Hermite polynomials to quantify the effects of variations in the load resistance, the eccentricity (distance between the center of mass and the elastic axis), and the nonlinear coefficients of the springs on the harvested power and the pitch and plunge amplitudes.The results show that the harvested power can be mostly affected by the eccentricity.Moreover, decreasing the nonlinear coefficient of the torsional spring results in a decrease in the pitch amplitude and an increase in the plunge amplitude and hence the harvested power.These results give guidance for optimizing and assessing the uncertainty in the performance of piezoaeroelastic energy harvesters.

Fig. 4 .
Fig.4.Probability density functions of the harvested power when using the Monte Carlo simulations and the first, second, fourth, seventh, and eighth non-intrusive PCE orders.

Fig. 6 .
Fig.6.Probability density functions of the pitch amplitude when using the Monte Carlo simulations and the first, second, fourth, seventh, and eighth non-intrusive PCE orders.
Probability density functions of the plunge amplitude when using the Monte Carlo simulations and the first, second, fourth, seventh, and eighth non-intrusive PCE orders.

Table 2 .
Sensitivity analysis of the harvested power and the plunge and pitch amplitudes to global variations in the load resistance, the eccentricity, and the nonlinear coefficients of the springs