Preliminary Study of Optimum Piezoelectric Cross-Ply Composites for Energy Harvesting

Energy harvesting devices based on a piezoelectric material attached to asymmetric bistable laminate plates have been shown to exhibit high levels of power extraction over a wide range of frequencies. This paper optimizes for the design of bistable composites combined with piezoelectrics for energy harvesting applications. The electrical energy generated during state-change, or “snapthrough,” is maximized through variation in ply thicknesses and rectangular laminate edge lengths. The design is constrained by a bistability constraint and limits on both the magnitude of deflection and the force required for the reversible actuation. Optimum solutions are obtained for differing numbers of plies and the numerical investigation results are discussed.


Introduction
Energy harvesting which converts ambient mechanical vibrations into electrical energy is an area of considerable research interest and has received extensive attention in the past decade.A variety of methods have been considered including inductive, capacitive, and piezoelectric materials [1][2][3].In many cases harvesting devices have been designed to operate at resonance to optimize the power generation, for example, simple linear cantilever beam configurations.However, ambient vibrations generally exhibit multiple time-dependent frequencies which can include components at relatively low frequencies.This can make typical linear systems inefficient or unsuitable; particularly if the resonant frequency of the device is higher than the frequency range of the vibrations it is attempting to harvest.In order to improve the efficiency of vibrational energy harvesters, recent work has focused on exploiting nonlinearity for broadband energy harvesting.Encouraging results [2] have been obtained using nonlinear or bistable cantilevered beams.Stanton et al. [2] modeled and experimentally validated a non-linear energy harvester using a piezoelectric cantilever.An end magnet on the oscillating cantilever interacts with oppositely poled stationary magnets, which induces softening or hardening into the system and allows the resonance frequency to be tuned.This technique was shown to outperform linear systems when excited by varying frequencies.However, such a system would require an obtrusive arrangement of external magnets and could generate unwanted electromagnetic fields.An alternative method has been recently found where a piezoelectric element is attached to bistable laminate plates with 2n plies and a total (T) layup of [0 n /90 n ] T to induce large amplitude oscillations [3].Such harvesting structures have been shown to exhibit high levels of power extraction over a wide range of frequencies.This arrangement can be designed to occupy a smaller space and is potentially more convenient and portable for broadband energy harvesting.
Bistable composites have been extensively studied for the development of morphing or adaptive structure concepts [4][5][6][7].When a composite laminate has an asymmetric stacking sequence the resulting mismatch in thermal expansion coefficients between plies leads to a thermally induced strain.This leads to the laminate developing a curved deformation as it is cooled from its high temperature cure cycle.Under certain geometric conditions the thermal strain can lead to the development of two stable equilibrium states ("bistability").Such structures are of interest for shape-change applications since the "snap-through" between stable-states results in a large deflection that does not require continuous energy input to be maintained.Figure 1 shows an example of this behavior for a square [0 n /90 n ] T laminate.For a low ratio of edge length to thickness only a saddle-shaped single stable state is observed, point A, with x-and y-curvatures of equal magnitude in opposite out-of-plane directions.As the ratio increases the solution bifurcates, point B. Beyond this point, two approximately cylindrical stable states are observed, points C and D, while the saddle state becomes unstable (dashed line).When a flexible piezoelectric material is attached to the laminate surface and the structure is repeatedly actuated between the stable states, the large shape changes have the potential to generate electrical energy by the direct piezoelectric effect.
Betts et al. [8] presented optimization of bistable laminates by deriving an analytical solution to modeling [0 n /90 n ] T as well as non-[0 n /90 n ] T layups.The primary focus of the paper was to develop shape changing structures whereby a piezoelectric strain was used to induce snapthrough of the bistable laminate between stable states 1 and 2. By variation of ply orientations and using a nonuniform laminate geometry, laminates were optimized for differing stiffness characteristics whilst satisfying bistability and deflection constraints.This work was extended to include piezoelectric layers into the existing analytical model for the design of laminates constrained by actuation voltage limits [9].This paper adapts this work and introduces optimization of bistable laminates combined with piezoelectric material for energy harvesting applications.The electrical energy generated during state-change (e.g., C to D in Figure 1) is maximized through variation in ply thicknesses and rectangular laminate edge lengths.A bistability constraint and limits on both the magnitude of deflection and the force required for the reversible actuation are specified to find a practical and achievable solution.Optimum solutions are obtained for differing numbers of plies and the numerical investigation results are discussed.

Design Problem
2.1.Problem Formulation.In this section we introduce the optimization formulation for the design of a bistable composite laminate for maximum piezoelectric energy harvesting capability as follows.
Maximize: Electrical energy generated by two orthogonal piezoelectric layers.
Subject to: The piezoelectric laminate must be bistable.
The stroke of actuation, w def , must be less than a maximum value.
The force F, required for state change from state 1 to state 2, and the reverse change from state 2 to state 1 must be less than a maximum specified value.
Variables: Piezoelectric layer and laminate ply orientations, θ p and θ respectively.

Nonuniform ply thicknesses.
Aspect ratio, L x /L y , of the laminate.

Laminate Configuration.
The laminate geometry has an orthogonal layup adhering to the following set of design rules as illustrated in Figure 2: (1) two piezoelectric layers, poled in the longitudinal direction (referred to as the d 33 direction), on the top and bottom surfaces of the laminate positioned 90 • apart; (2) an even number of plies, with pairs of plies about the laminate midplane positioned 90 • apart; (3) equal thicknesses, t p , for both piezoelectric layers; (4) equal ply thicknesses about the laminate midplane, t 1 , t 2 , . . ., t n ; (5) nonuniform edge lengths, L x and L y , common to both the laminate and piezoelectric layers; (6) each ply of the laminate is made from the same material.
The laminate geometry above is selected for three reasons.Firstly the piezoelectric layers are positioned orthogonally to align with the major curvatures of the two stable  shapes, maximizing the stress in the piezoelectric elements and, hence, the mechanical to electrical energy conversion.
Secondly the laminate stacking sequence allows maximizing the useful deflection between states while providing scope for tailoring the directional properties.Finally, the restricted orthogonal stacking sequence reduces the design search space by excluding general asymmetric layups which do not exhibit the high deflection bistable characteristics.It is noted that the laminate (excluding the piezoelectric layers) can be defined in terms of plies above the midplane, referred to here as the half-laminate, (θ 1 and θ 2 in Figure 2).We initially consider a rectangular cross-ply design, [0 n /90 n ] T , as studied experimentally by Arrieta et al. [3] This laminate exhibits equal out-of-plane displacement at all four corners and is actuated in the configuration shown in Figure 3.The laminate is held at all four corners (zero z-displacement).A mechanical actuator is attached to the centre point of the laminate surface such that displacement applied in the z-direction (out-of-plane) induces a statechange, simulating the ambient mechanical motion.The out-of-plane displacement, referred to here as the "stroke of actuation," is the change in maximum out-of-plane deflection from one unloaded stable state to the second, w def , measured at the geometric centre of the laminate.

Objective Function.
The objective of this study is to maximize the electrical energy generated by the mechanical deflection of a bistable piezoelectric-composite laminate.We model this based on the static states of the system.The objective function is thus the electrical energy, U, induced by stress in the piezoelectric material induced in each stable state: where C is the capacitance and V is the open circuit output voltage, as defined by the following equations: where d i j is the effective piezoelectric strain constant (charge per unit force), g i j is the effective piezoelectric voltage constant (electric field per unit stress), σ is the stress, A is the surface area of the piezoelectrics, and t p is the thickness of the layer.Substituting ( 2) into (1) gives: where the terms have been split into the fixed material properties of the piezoelectric (d i j and g i j ), the stress in the material which is a function of the laminate curvature, and the material geometry (A and t p ).When attached to the surface of a curved laminate the stress in the piezoelectric material varies across the surface area as a function of the strains induced by the mismatched thermal properties.The piezoelectric material, poled in the longitudinal direction in this work, is strained in both the longitudinal direction (referred to as the 33 direction) and the transverse direction (31 direction).We therefore consider the contributions of both components.The stresses in the x-and y-directions are defined by the following equations, where the Q's are transformed stiffness terms dependent on the material orientations and ε's are the in-plane strains.These values vary across the x-y plane and are therefore integrated over the piezoelectric surface area and substituted into (3).We then consider the total electrical energy induced by a single stable state to be the sum of the longitudinal and transverse components of both the top and bottom piezoelectric layers.Each of these terms uses the associated 33 and 31 material properties as illustrated in Figure 4.The total electrical energy, defining the objective function is therefore, the sum of the four terms associated with state 1, and the four terms associated with state 2, (5): (5)

Constraints on the Design
2.4.1.Bistability.The asymmetry of laminates does not guarantee bistability.When a laminate geometry has an edge length to thickness ratio below a critical value there exists only a single saddle solution [10].This critical ratio is a function of the stacking sequence, temperature, material properties, and aspect ratio.Furthermore, the critical value is increased by the presence of the piezoelectric layers due to their added stiffness.We introduce a constraint to ensure bistability.Noting that the x-and y-curvatures of the monostable saddle shape have equal magnitude in opposite out-of-plane directions, the sum of these curvatures is zero.Beyond the bifurcation point, see Figure 1, the sum of the curvatures is either positive or negative depending on the associated cylindrical state.We therefore impose the constraint a + b > 0, where a and b are the x-and y-curvature coefficients, respectively, such that only bistable composites are considered feasible: Close to the bifurcation point the solution is highly sensitive to manufacturing imperfections and uncertainties in material data [11].This issue can easily be addressed by adding a small constant to the constraint, that is, a + b > X, where X is a "distance" away from the bifurcation point.Our previous studies have shown that this constraint could also be used as a minimum out-of-plane deflection constraint [8].

Stroke of Mechanical Actuation.
The stroke of actuation refers to the distance travelled by the mechanical actuator in the z-direction, Figure 3 and a practical limit is imposed.This is modeled by fixing the displacement at the centre with the maximum deflection being at the corners of the laminate.Noting that the out-of-plane displacement, w of the curved laminate shape is approximated by quadratic polynomial of (7): the corner deflection between states, w def , at x = L x /2, y = L y /2, can be expressed as (8) where the subscript denotes the associated stable state.
Due to the orthogonality of the laminate the deflection between states is the same at all four corners.Furthermore, the curvatures a, b, and c can be expressed in terms of the first state alone due to the two states having curvatures of equal magnitude.The constraint on the stroke, less than w max , can therefore be reduced to (9): Actuation force constraint:

Analytical Model
3.1.Unloaded Laminate Shapes.The analytical model used to calculate the unloaded shapes of asymmetric laminates of arbitrary layup was introduced by Dano and Hyer [5] and is a nonlinear extension to classical laminated plate theory.The co-ordinate system used is that defined in Figure 3, where the geometric centre of the laminate sits at the origin and the ply orientations are measured from the x-axis.The out-ofplane displacement in the z-direction, w is assumed to be of the form of (7).The midplane strains, including geometrical nonlinearity according to the von Karman hypothesis, are defined as where u 0 and v 0 are the in-plane displacements in the x-and y-directions, respectively.The midplane strains are approximated by the third order polynomials.Dano and Hyer [5] found that terms with powers of x and y that sum to an odd number are always zero.Therefore the form of the midplane strains can be reduced to the polynomials of ( 12): Using ( 7), (11), and ( 12) and introducing the additional shape coefficients e 9-11 resulting from integration of the midplane strains, expressions for the in-plane displacements u 0 and v 0 can be determined: (13) The total strain energy of the laminate, W, can then be expressed as the integral of strain energy density over the volume of the laminate.
where c i jkl is elastic constant, α i j is constant relating to the thermal expansion coefficients, L x and L y are the planform side lengths of the laminate, H is the total laminate thickness, ΔT is the temperature change from the cure temperature, and ε i j 's and ε kl 's are the total strains defined as: where z is the out-of-plane distance from the laminate midplane and a, b, and c are out-of-plane shape coefficients defined by (7).Expansion of ( 14) results in an expression for the total energy which is a function of the material and geometric properties, the temperature change from the laminate cure and the set of shape coefficients a, b, c, e 1 . . .e 11 .
For equilibrium, the minimum energy states require: where the k i 's are the set of shape coefficients.For the laminate stacking sequence of Figure 2, with square laminate edge lengths, Betts et al. [9] derived an analytical solution to (16) with unloaded orthogonal piezoelectric layers.This provides a fast and reliable method of calculating all stable shapes of asymmetric square laminates, necessary for optimization.In this work we use this analytical solution for a square laminate of edge lengths L x as the initial solution for rectangular laminate shapes and the Newton-Raphson method is employed to solve the system of (16).

Actuation.
In the arrangement shown in Figure 3, as used by Arrieta et al. [3], the out-of-plane force F is applied at the centre of the laminate with the four corners held.We model this by fixing the origin and the maximum deflection is at the corners where one quarter of the force is applied at each corner, Figure 5.
The out-of-plane deflection, δw is evaluated using ( 7) where w a denotes unloaded out-of-plane deflection and w b the deformed out-of-plane deflection: The work done in deflecting the laminate, W F , is hence The total energy, W tot , ( 19) is minimized, where W represents the strain energy of (14): Equation ( 20) is solved using a simple iterative scheme.The exact unloaded laminate shape is obtained first as described in Section 3.1 and the force is incrementally increased until the state transition.Figure 6 exhibits an example of determining a snap-through force where the applied force and the corresponding curvatures are plotted to show the state transition.Starting with the unloaded laminate at A there is a large positive x-curvature for state 1. Between points A and B the force is increased.From point B the curvature is seen to change dramatically with x-curvature becoming very small and the y-curvature becoming large and negative at point C.This point marks the transition from state 1 to state 2. As the load is then incrementally removed the laminate curvatures tend towards point D which represents the unloaded second state.This process can be repeated in the opposite out-of-plane direction to reset the laminate and calculate the second snap force.These two values define the force constraints for optimization.

Results
This section presents the results of a design study.The optimization problem is solved using Matlab's sequential quadratic programming, fmincon.All examples use M21/T800 material properties [12] for the laminate and M8557-P1 MFC properties for the piezoelectric layers [13].

Influence of Aspect Ratio for a Single Laminate.
Prior to optimization, we first investigate the influence of aspect ratio on the strain energy components of (5) using [0/0/90/90] T laminates of aspect ratio 1 and 5 (Figure 7).A breakdown of the objective function values is shown in Table 1.For the laminate with an aspect ratio of 1 (Figure 7(a)), two equal and opposite shapes are observed and each state contributes 50% of the total objective function.Of each half, the vast majority is contributed by the longitudinal direction of the piezoelectric aligned with the major curvature, U t x for state 1 and U b y for state 2. Interestingly the same pattern is observed for the aspect ratio of 5 where the two stable states have very different displacements, Figure 7(b).Despite the appearance, states 1 and 2 exhibit the same magnitude of curvature in opposite out-of-plane directions, thus the vast majority of the electrical energy is provided by the longitudinal direction of the piezoelectric layer aligned with the major curvature.Our detailed investigation revealed that this is consistently observed for all cross-ply laminates.

Optimization of Cross-Ply
Laminates.This section presents optimization of a [0/0/90/90] T laminate to determine ply thickness, t, and the laminate aspect ratio, L x /L y .The piezoelectric layers are oriented at 0 • and 90 • and have a fixed thickness of 0.2 mm.The surface area of the laminate and the piezoelectric layers is fixed to 0.04 m 2 .The full design space for this problem constrained by bistability, a maximum deflection of 25 mm and a maximum reversible snap force of 5 N is shown in Figure 8.The colored contours represent the objective function value bounded by each of the constraints.The outer limits of the design space (the white "monostable" area) are due to the bistability constraint and are observed to be at a maximum half-laminate thickness of 0.866 mm (with an aspect ratio 1), and a maximum aspect ratio of 6.85 (at a half-laminate thickness of 0.196 mm), Figure 8.The left hand boundary of the monostable region defines the minimum laminate thickness to achieve bistability.As the thickness of the piezoelectric is fixed, the thickness ratio of the laminates to piezoelectrics becomes too small beyond this minimum limit and the stiffness contribution of the piezoelectrics to the piezo-composite laminate becomes significant, leading to the loss of bistability.The right hand boundary represents the bistability bifurcation limit due to low edge to thickness ratio.As shown in Figure 1, there is only one saddle equilibrium state for low edge to thickness laminates.For the optimization study presented here, the laminate area is fixed which limits the maximum achievable edge length.As the laminate thickness increases on the x-axis in Figure 8, the edge to thickness ratio decreases and reaches the right hand boundary beyond which the laminates can only have one equilibrium state.There are three grey areas marking the infeasible regions.The lightest region indicates laminates which violates the actuation stroke constraint, that is, they require a deflection greater than the allowable maximum.This constraint violation, however, does not influence the optimum designs as the objective function values of these solutions are relatively low.
The other two constraints, snap-through forces F 1 and F 2 , required to induce state transition, however, are critical and determine the optimum solutions.These are two darker grey regions in Figure 8 where the objective function values are maximum.The laminates in these regions require a greater force to induce the state change than the available force.For the case of a [0/0/90/90] T laminate, the snap force from state 1 to 2 is completely enclosed by the snap force from state 2 to 1.This is a result of the aspect ratio where the displacements between two states can vary significantly for high aspect ratio (e.g., Figure 7).The difference between the snap force requirements is more prominent for low thickness laminates in Figure 8 and indeed, these regions meet at aspect ratio 1 indicating that the force requirement from states 1 to 2 and 2 to 1 are identical for square laminates.The lower limit of the constraint boundary indicates that the snap force requirement is reduced as the thickness of the laminates decreases.As the thickness increases, the increased stiffness of the laminates reduces the curvature thus requiring less snap force and displacement.This gives rise to the upper constraint boundary.
Two optimal solutions are obtained.The first, indicated as a white circle in Figure 8 is a local optimum with an objective function value of 0.012 J.This laminate has an aspect ratio of 1 and half-laminate thickness of 0.245 mm.The second, the black circle, has an aspect ratio of 1 and halflaminate thickness of 0.731 mm, giving an objective function value of 0.0295 J.The second solution is the global optimum.

Discrete Ply Thicknesses.
The problem so far assumed that the ply thickness variation is continuous.However in practice, the number of plies and, thus, the laminate thickness are discrete.For a typical single ply thickness of 0.125 mm, we can reconsider the design space of Figure 8 as a [0 n /90 n ] T laminate with n ranging from 1 to 6, that is, 2to 12-ply laminates.With this added constraint the discrete optimum solution for each number of plies are marked by white squares in Figure 8 and summarized in Table 2.The solutions tend to be for laminates with aspect ratio 1 unless the square laminates violate a constraint.The active constraint in this case was the mechanical force limit.The global solution was found to be a square laminate with n = 6, whose objective function is an order of magnitude greater than of the n = 1 solution.

Conclusions
This paper presented a modeling and optimization formulation of bistable piezo-composite laminates for energy harvesting applications.The method allows design param- eters such as harvester aspect ratio and device thickness to be optimized based on realistic design constraints such as snap-through forces, displacements, and the influence of the stiffness of the piezoelectric elements on bistability of the laminate.This preliminary investigation of optimum crossply laminates has found that despite the high deflection of high aspect ratio laminates, the square laminates are the optimum with highest electrical energy.The mechanical force required to induce state change was found to be the critical factor in determining the optimum solution whilst the deflection constraint was observed to be inactive.The discrete ply thickness and its influence on the solutions were also examined.Based on this preliminary study, continuing research is extending the optimization methodology to include additional design constraints such as the maximum stress in the piezoelectric material to prevent depolarization.Further design parameters are also considered such as choice of piezoelectric material, piezoelectric area, and configurations and stacking sequence of non-cross-ply laminates.

Figure 3 :
Figure 3: Arrangement for actuation of a rectangular cross-ply laminate.

1 2Figure 4 :
Figure 4: Longitudinal and transverse components of the total electrical energy, U. Note: subscript denotes the associated direction; superscript denotes the top (t) and bottom (b) piezoelectric layer.

Figure 6 :
Figure 6: Changing curvature and snap-through behavior of a [0 n /90 n ] T laminate with applied force.

Table 2 :
Optimum solution for discrete laminate thicknesses.