1 INTRODUCTION
Beam formulations are widely used for the numerical analysis of one dimensional piezoelectric structure (^{Marinkovic and Marinkovic, 2012}). Analytical formulations (^{Crawley and de Luis, 1987}; ^{Abramovich and Pletner, 1997}; ^{Crawley and Anderson, 1990}) and finite elements (^{Bendary et al., 2010}; ^{Kumar and Narayanan, 2008}; ^{Sulbhewar and Raveendranath, 2014a}) based on the Euler-Bernoulli beam theory can be effectively used for the analysis of thin and slender piezoelectric smart beams. However, Euler-Bernoulli theory neglects the deformation due to shear and hence not suitable for thick and short beams. Sandwich Beam Theory (SBT) based analytical formulation (^{Zhang and Sun, 1996}) and finite element formulations (^{Benjeddou et al., 1997}, ^{2000}; ^{Raja et al., 2002}) considered the thicker core as a Timoshenko beam and the relatively thinner faces as Euler-Bernoulli beams. However, SBT is not suitable for short beams with thick piezoelectric layers. Timoshenko beam formulations based on the First-order Shear Deformation Theory (FSDT) consider constant shear strain across the beam cross-section. Analytical formulations (^{Abramovich, 1998}; ^{Aldraihem and Khdeir, 2000}; ^{Khdeir and Aldraihem, 2001}) and finite elements (^{Robbins and Reddy, 1991}; ^{Shen, 1995}; ^{Narayanan and Balamurugan, 2003}; ^{Ray and Mallik, 2004}; ^{Neto et al., 2009}) based on FSDT are widely used in the literature for the analysis of piezoelectric smart structures.
Accuracy of the conventional FSDT-based piezoelectric beam finite elements (^{Shen, 1995}; ^{Narayanan and Balamurugan, 2003}; ^{Ray and Malik, 2004}; ^{Neto et al., 2009}) is adversely affected by the induced potential effects. These elements consider linear through-thickness distribution of electric potential which is actually nonlinear by virtue of the induced potential. The accuracy can be improved using assumed higher-order approximation of through-thickness electric potential (^{Jiang and Li, 2007}; ^{Kapuria and Hagedorn, 2007}; ^{Wang et al., 2007}; ^{Beheshti-Aval and Lezgy-Nazargah, 2012}, ^{2013}). However, assumed higher-order potential distribution in the formulation introduces additional nodal electrical degrees of freedom in the transverse direction and hence increases the computational cost. An alternate efficient way to include the higher-order induced potential in FSDT-based formulation is to use the consistent through-thickness potential distribution derived from the electrostatic equilibrium equation (^{Sulbhewar and Raveendranath, 2014b}).
Also, convergence of the conventional two-noded isoparametric FSDT-based piezoelectric beam element (^{Narayanan and Balamurugan, 2003}) depends on the extent of the extension-bending and bending-shear couplings. Recently, (^{Sulbhewar and Raveendranath, 2015}) proposed a novel FSDT piezoelectric beam finite element based on coupled polynomials for field variables which showed improved convergence. However, this element is not consistently accurate as the governing equations used to define coupled shape functions in this formulation are based on the assumed linear through-thickness potential. The associate errors are prominent for beams with piezoelectrically dominant cross-sections and/or with higher piezoelectric coefficients. An ideal FSDT-based formulation which is accurate and efficient over all geometric and material configurations of the piezoelectric beam should incorporate the induced potential coupling along with other mechanical couplings at the field interpolation level itself.
In the present work, an attempt is made to develop a novel FSDT piezoelectric beam formulation which is consistently accurate and efficient throughout the applicable range of geometric and material properties. The governing equations are derived using the variational formulation based on FSDT in conjunction with the consistent through-thickness potential. The relations established between field variables are used to define coupled quadratic polynomials for axial displacement (u_{o} ) and section rotation (θ), having contributions from the assumed cubic polynomial for transverse displacement (w_{o} ) and assumed linear polynomials for layerwise electric potential variables . The shape functions based on these polynomials efficiently handle change in stiffness due to the induced potential along with bending-extension and bending-shear couplings, in an efficient manner. Comparison of results from the present and the conventional formulations against ANSYS 2D benchmark simulation results proves the improved accuracy of the present formulation over the conventional formulations. Convergence studies are carried out to prove the improved convergence characteristics of the present FSDT element over the conventional isoparametric FSDT beam elements. It is noteworthy that owing to the fully coupled polynomial representation for section rotation and coupled quadratic term in the interpolation polynomial for axial displacement and transverse electric potential, the improved performance has been achieved with the same number of nodal degrees of freedom as of conventional two-noded isoparametric FSDT-based piezoelectric beam element.
2 THEORETICAL FORMULATION
An equivalent single layer (ESL) FSDT model for mechanical fields and a layerwise model for electric potential (φ) are employed for the proposed formulation. Consider a general multilayered extension mode piezoelectric smart beam with total number of layers n, as shown in Figure 1. The layers can be host layer(s) of conventional material or bonded/embedded layers of piezoelectric material. The beam layers are assumed to be made up of isotropic or specially orthotropic materials with perfect bonding among them. Top and bottom faces of piezoelectric layers are fully covered with electrodes. Mechanical and electrical quantities are assumed to be small enough to apply linear theories of elasticity and piezoelectricity and assumptions of beam theory apply.
2.1 Reduced Constitutive Relations
For a general piezoelectric smart structure, the elastic Cij (i, j = 1...6), piezoelectric ekj(k = 1,2,3, j=1....6) and dielectric ∈k (k = 1,2,3) constants relate the mechanical and electrical variables through the three-dimensional constitutive relations. An extension mode piezoelectric smart beam with axes of material symmetry parallel to the beam axes is considered here. For extension mode, the transversely poled piezoelectric material is subjected to the transverse electric field. For such a beam, the general constitutive relations are reduced according to beam theory, which are given as (^{Sulbhewar and Raveendranath 2014 b}):
where (i=1....number of piezoelectric layers), (k=1.....n). σ, τ, ε, γ, D and E denote the axial stress (N/m2), shear stress (N/m2), normal strain, shear strain, electric displacement (C/m2) and electric field (V/m), respectively. The constants (i = 1,5), ẽ and denote reduced elastic (N/m ^{2)}, piezoelectric (C/m ^{2)} and dielectric (F/m) properties, respectively.
2.2 Mechanical Displacements and Strains
The mechanical displacement fields in the longitudinal and transverse directions for FSDT are given as (^{Narayanan and Balamurugan, 2003}):
u _{0}(x) and w _{0}(x) are the centroidal axial and transverse displacements, respectively. θ is the section rotation of the beam. Dimensions L,b,h denote the length, width and the total thickness of the beam, respectively.
Axial and shear strain fields are derived using usual strain-displacement relations as:
where ()' denotes derivative with respect to x.
2.3 Electric Potential and Electric Field
The layerwise two dimensional electric potential φ _{i} (x,z) takes the values of Φ _{i + 1} (x) and Φ _{i} (x) at the top and bottom faces of ^{ ith } piezolayer, respectively as shown in Figure 1. The through-thickness distribution of electric potential φ _{i} (x,z) consistent with FSDT is used (^{Sulbhewar and Raveendranath, 2014b}):
where
The first two terms of expression ^{(6)} describe the conventional linear part in which and are the mean and difference, respectively, of the top and bottom surface potentials of the ^{ ith } piezoelectric layer. The quadratic term represents the coupling between curvature strain and electric potential which constitutes induced potential.
The layerwise electric field ( ^{ Ei } _{z} ) is obtained from equation ^{(6)} as (^{Sulbhewar and Raveendranath, 2014b}):
3 VARIATIONAL FORMULATION
The formulation is based on Hamilton's principle which implicitly takes care of natural boundary conditions. It is expressed as (^{Chee et al., 1999}):
where, K=kinetic energy, H=electric enthalpy density function for piezoelectric material and mechanical strain energy for the linear elastic material and W=external work done.
3.1 Variation of Electromechanical/Strain Energy
The electromechanical/strain energy variation of the piezoelectric smart beam is given as (^{Chee et al., 1999}):
Substituting values of axial strain (ε _{x} ), shear strain (γ _{xz} ), electric field ( ^{ Ei } _{z} ) from equations ^{(4)}, ^{(5)}, ^{(7)} and using them along with constitutive relations given by equation ^{(1)} in expression ^{(9)}; the total variation on the potential energy of the smart beam is given as:
3.2 Variation of Kinetic Energy
Total kinetic energy of the beam is given as (^{Chee et al., 1999}):
where ρ _{k} is the mass density of ^{ kth } layer in kgm ^{-3} and (k =1...n). Substituting values of u and w from equations ^{(2)} and ^{(3)} and applying variation, to derive at:
where ^{(˙)} denotes ∂/∂t.
3.3 Variation of Work of External Forces
Total virtual work of the structure can be defined as the product of virtual displacements with forces for the mechanical work and the product of the virtual electric potential with the charges for the electrical work. The variation of total work done by external mechanical and electrical loading is given by (^{Chee et al., 1999}):
in which ^{ fV, fS, fC } are volume, surface and point forces, respectively. q _{0} and S _{φ} are the charge density and area on which charge is applied.
4 DERIVATION OF COUPLED FIELD RELATIONS
The relationship between the field variables is established here using static governing equations. For static conditions without any external loading, the variational principle given in equation ^{(8)} reduces to (^{Sulbhewar and Raveendranath, 2015}):
Applying variation to the basic variables in equation ^{(10)}, the static governing equations are obtained as:
where
Assuming that the higher order continuous derivatives of variables appearing in the governing equation ^{(17)} exist, we can write:
Using equations ^{(15)} and ^{(18)}, we can write the relationship of axial displacement (u _{0}) with transverse displacement (w _{0}) and electric potential variable as:
where β_{1} = ^{ Au } _{2}/ ^{ Au } _{1} and β ^{i} _{2} = ^{ Aui } _{3}/ ^{ Au } _{1}.
From equations ^{(16)}-^{(19)}, we can write the relationship of section rotation (θ) with transverse displacement (w _{0}) and electric potential variable as:
These equations for u"0 and θ are used in the next Section to derive coupled polynomial expressions for the field variables. It is clear that the coupling coefficients βj(j = 1,2, 3, 4) which depend on geometric and material properties of the beam, relate all the field variables by properly accommodating bending-extension, bending-shear and induced potential couplings. It is noteworthy that the constants Anm(m = 1;2 3; 4) (n = u, θ w) appearing in equations ^{(15)}-^{(17)}, which are used to define the coupling coefficients βj(j = 1,2, 3, 4) are different from those given in the (^{Sulbhewar and Raveendranath, 2015}). The constants Anm in the present formulation contain additional stiffness terms (shown in curly braces) due to the induced potential effects. It may be noted that this induced stiffness is proportional to which bears the same unit ( ^{ N/m2 } ) as of elastic modulus . Hence, the quantity may be termed as 'induced modulus'.
5 FINITE ELEMENT FORMULATION
Using the variational formulation described above, a finite element model is developed here. The two-noded beam element considered here is based on FSDT with layerwise electric potential in the transverse direction. There are three mechanical variables in the formulation namely, u _{0}, w _{0}and θ and layerwise electric potential variables where (i =1.....number of piezoelectric layers in the beam).
The equations ^{(19)} and ^{(20)}, derived using the governing equilibrium equations, demand continuous third order derivative of w _{0} and first order derivative of . Hence, in terms of the natural coordinate ξ, a cubic polynomial for transverse displacement w0 and linear polynomials for layerwise electric potential variable are assumed as given in equations ^{(21a)} and ^{(21b)}, respectively. The transformation between the local coordinate ξ and the global coordinate x along the length of the beam is given as ξ = [2(x - x _{1}) / (x _{2} - x _{1})] - 1 and (x _{2} - x _{1}) - l, length of the beam element.
Using equations ^{(21a)} and ^{(21b)} in equation ^{(19)} and integrating with respect to ξ, we get the coupled polynomial for axial displacement u _{0} as:
It is noted that the coupled quadratic term in equation ^{(22)} contains contributions from w _{0} and fields and does not bring in any additional generalized degree of freedom.
Substituting equations ^{(21a)} and ^{(21b)} in equation ^{(20)}, the coupled polynomial expression for the section rotation θ is derived as:
Equation ^{(23)} interpolates θ by purely coupled terms with contributions from w _{0} and fields.
It is noteworthy that equations ^{(22)} and ^{(23)} take care of extension-bending, bending-shear and induced potential couplings in a variationally consistent manner with the help of coupled terms present in the description of axial displacement and section rotation.
Using equations ^{(21)}-^{(23)}, the coupled shape functions [ ^{ Nu } _{m} ](m = 1..8), [ ^{ Nw } _{m} ](m = 1..6), [N ^{θ} _{m} ](m = 1..6) and [N ^{φ} ^{i} _{m} ](m = 1..2) which relate the field variables to their nodal values as given in equation ^{(24)} are derived by usual method. The expressions for shape functions are given in ^{Appendix}.
As noted from the equation ^{(24)}, while employing quadratic polynomials for axial displacement u _{0} and section rotation θ in the present FSDT formulation, the number and type of nodal variables are maintained the same as of the conventional isoparametric FSDT formulation.
The variation on basic mechanical and electrical variables can now be transferred to nodal degrees of freedom. Substituting equation ^{(24)} in equations ^{(10)}, ^{(12)}, ^{(13)} and using them in equation ^{(8)}, the following discretized form of the model is obtained:
where M is mass matrix, K_{uu}, K_{uφ}, K_{φu}, K _{φφ} are global stiffness sub-matrices. U, Φ are the global nodal mechanical displacement and electric potential degrees of freedom vectors, respectively. F and Q are global nodal mechanical and electrical force vectors, respectively. The matrix equations are now solved according to electrical conditions (open/closed circuit), configuration (actuator/sensor) and type of analysis (static/dynamic).
6 NUMERICAL EXAMPLES AND DISCUSSIONS
The software implementation of the present formulation has been carried out in MATLAB environment. The accuracy and efficiency of the proposed FSDT finite element are tested here for static (actuation/sensing) and modal (open/closed circuit) analyses and its performance is compared against the conventional two-noded isoparametric FSDT piezoelectric beam finite element available in the literature. The following designations are used:
FSDT-Coupled: The present formulation which uses the coupled polynomials (cubic for w _{0} given by equation ^{(21a)}, coupled quadratic for u _{0} given by equation ^{(22)}, coupled quadratic for θ given by equation ^{(23)} and linear for given by equation ^{(21b)}) for interpolation of field variables and layerwise consistent through-thickness potential (coupled quadratic approximation in z direction given by equation ^{(6)}).
FSDT: The conventional FSDT formulation of (^{Narayanan and Balamurugan, 2003}) which uses independent polynomials for field interpolation (linear for u _{0,} w _{0,} θ and ) and layerwise assumed linear through-thickness potential.
ANSYS 2D: For a comparative evaluation of the above FSDT formulations, benchmark solutions have been obtained from a refined two-dimensional analysis using ANSYS finite element software, for which PLANE 183 elements are used to mesh conventional material layers, while PLANE 223 elements are used to mesh piezoelectric material layers.
6.1 Example 1: A Symmetric Bimorph Beam
A bimorph cantilever beam with oppositely poled piezoelectric layers as shown in Figure 2 is considered here. In order to study the effect of material properties on the performance of the FSDT elements, the following materials are used while the geometry is fixed (h = 10mm, L = 100mm).
PVDF (^{Sun and Huang, 2000}):
E = 2Gpa, υ = 0.29, e - 0.046Cm-2, ∈ = 0.1062 × 10-9Fm-1, ρ = 1880kgm-3
PZT 2 (^{EFunda.com, 2014})
{C11, C12, C13, C22, C23, C33, C44, C55, C66} = {134.87, 67.89, 68.09, 134.87, 68.09, 113.30, 22.22, 22.22, 33.44}GPa
{e31, e32, e33} = {-1.8160, -1.8160, 90506}Cm-2 {∈1, ∈2, ∈3} = {8.7655, 8.7655, 3.9843} × 10-9Fm-1
ρ = 7600kgm-3
PZT 4 (^{EFunda.com, 2014})
{C11, C12, C13, C22, C23, C33, C44, C55, C66} = {139, 77.84, 74.28, 139, 74.28, 115.41, 25.64, 25.64, 30.58}GPa
{e31, e32, e33} = {-5.2028, -5.2028, 15.0804}Cm-2 {∈1, ∈2, ∈3} = {1.3060, 1.3060, 1.1510} × 10-8Fm-1
ρ = 7500kgm-3
PZT-5H (^{Kapuria and Hagedorn, 2007}):
{C11, C12, C13, C22, C23, C33, C44, C55, C66} = {126, 79.5, 84.1, 126, 84.1, 117, 23, 23, 23.25}GPa
{e31, e32, e33} = {-6.5, -6.5, 23.3}Cm-2 {∈1, ∈2, ∈3} = {1.503, 1.503, 1.3} × 10-8Fm-1
ρ = 7500kgm-3
PZT 5A (^{EFunda.com, 2014})
{C11, C12, C13, C22, C23, C33, C44, C55, C66} = {120.35, 75.18, 75.09, 120.35, 75.09, 110.87, 21.05, 21.05, 22.57}GPa
{e31, e32, e33} = {-5.3512, -5.3512, 15.7835}Cm-2 {∈1, ∈2, ∈3} = {1.5317, 1.5317, 1.5052} × 10-8Fm-1
ρ = 7550kgm-3
PZT 8 (^{EFunda.com, 2014})
{C11, C12, C13, C22, C23, C33, C44, C55, C66} = {146.88, 81.09, 81.05, 146.88, 81.05, 131.71, 31.35, 31.35, 32.89}GPa
{e31, e32, e33} = {-3.8754, -3.8754, 13.9108}Cm-2 {∈1, ∈2, ∈3} = {1.1422, 1.1422, 0.8854} × 10-8Fm-1
ρ = 7600kgm-3
G1195N (^{Peng et al., 1998}): E = 63GPa, υ = 0.3, d31 = 254×10-12mV-1, e31 = E×d31Cm-2, ∈3 = 15×10-9Fm-1, ρ = 7600kgm-3
For a comparative evaluation of accuracy of various FSDT-based formulations, converged results from a refined mesh of 40 elements have been used. Converged results from an ANSYS 2D simulation with a mesh of 100×10 elements are used as a benchmark.
Static Analysis-Actuator Configuration: In this configuration, the interface of the bimorph is grounded and the voltages of ±10volts are applied on the free surfaces. Table 1 and ^{2} show the results for the tip deflection and the maximum axial stress developed in the bimorphs of different materials, respectively. Also, the associated absolute errors calculated with respect to ANSYS 2D benchmark solutions are presented in brackets. As seen from the tables, the conventional FSDT formulation (^{Narayanan and Balamurugan, 2003}) fails to produce consistently accurate results. The percentage errors increase with the modulus ratio (the ratio of induced modulus to elastic modulus) of the material. The present FSDT-Coupled element predicts accurate results for all the bimorphs, irrespective of the modulus ratio. This consistent performance of the present formulation can be attributed to the accommodation of induced potential effects through the coupled interpolation polynomials.
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
PVDF | 0.0100 | 0.0690 (0.291 %) | 0.0688 | 0.0688 (0.000 %) |
PZT-2 | 0.0589 | 0.1800 (1.408 %) | 0.1775 | 0.1774 (0.056 %) |
PZT-8 | 0.0687 | 0.2910 (1.677 %) | 0.2862 | 0.2861 (0.035 %) |
PZT-4 | 0.0773 | 0.3690 (1.906 %) | 0.3621 | 0.3620 (0.028 %) |
PZT-5A | 0.0849 | 0.5130 (2.090 %) | 0.5025 | 0.5023 (0.040 %) |
PZT-5H | 0.1751 | 0.8244 (4.328 %) | 0.7902 | 0.7899 (0.038 %) |
G1195N | 0.2710 | 0.7620 (6.723 %) | 0.7140 | 0.7137 (0.042 %) |
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
PVDF | 0.0100 | 0.0460 (0.648 %) | 0.0463 | 0.0463 (0.000 %) |
PZT-2 | 0.0589 | 5.1724 (4.174 %) | 5.3977 | 5.3977 (0.000 %) |
PZT-8 | 0.0687 | 8.4348 (4.810 %) | 8.8610 | 8.8618 (0.009 %) |
PZT-4 | 0.0773 | 10.000 (5.366 %) | 10.567 | 10.569 (0.019 %) |
PZT-5A | 0.0849 | 10.427 (5.868 %) | 11.077 | 11.077 (0.000 %) |
PZT-5H | 0.1751 | 16.492 (11.18 %) | 18.567 | 18.567 (0.000 %) |
G1195N | 0.2710 | 16.002 (15.99 %) | 19.048 | 19.048 (0.000 %) |
Static Analysis-Sensor Configuration: Here, the beam shown in Figure 2 is subjected to a tip load of 1000 N. The results for the tip deflection, potential developed at the mid-span and the maximum axial stress developed at the root of the bimorphs of different materials are tabulated in the Tables 3, ^{4} and ^{5}, respectively. The associated absolute errors (in percentage) with respect to ANSYS 2D benchmark solutions are presented in brackets. As seen from the tables, the present FSDT-Coupled consistently reproduces the ANSYS 2D simulation results for all the bimorphs, unlike the conventional formulation. The accuracy of the present FSDT-Coupled formulation is practically insensitive to the material properties of the beam.
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
PVDF | 0.0100 | 2000.0 (0.090 %) | 1998.2 | 2000.0 (0.090 %) |
PZT-2 | 0.0589 | 44.885 (1.456 %) | 44.241 | 44.268 (0.061 %) |
PZT-8 | 0.0687 | 44.065 (1.373 %) | 43.468 | 43.364 (0.239 %) |
PZT-4 | 0.0773 | 46.894 (1.661 %) | 46.128 | 46.061 (0.145 %) |
PZT-5A | 0.0849 | 62.147 (1.752 %) | 61.077 | 60.942 (0.221 %) |
PZT-5H | 0.1751 | 59.349 (3.757 %) | 57.200 | 57.156 (0.077 %) |
G1195N | 0.2710 | 53.180 (5.247 %) | 50.529 | 50.368 (0.319 %) |
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
PVDF | 0.0100 | 1612.2 (0.261 %) | 1608.0 | 1608.3 (0.019 %) |
PZT-2 | 0.0589 | 81.813 (1.411 %) | 80.675 | 80.675 (0.000 %) |
PZT-8 | 0.0687 | 58.055 (1.642 %) | 57.117 | 57.123 (0.011 %) |
PZT-4 | 0.0773 | 54.778 (1.835 %) | 53.791 | 53.795 (0.007 %) |
PZT-5A | 0.0849 | 57.410 (2.001 %) | 56.284 | 56.287 (0.005 %) |
PZT-5H | 0.1751 | 70.381 (3.871 %) | 67.758 | 67.758 (0.000 %) |
G1195N | 0.2710 | 105.55 (5.628 %) | 99.926 | 99.924 (0.002 %) |
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
PVDF | 0.0100 | 5.955 (0.750 %) | 6.000 | 5.993 (0.117 %) |
PZT-2 | 0.0589 | 5.886 (2.517 %) | 6.038 | 5.993 (0.745 %) |
PZT-8 | 0.0687 | 5.873 (2.861 %) | 6.046 | 5.994 (0.860 %) |
PZT-4 | 0.0773 | 5.862 (3.011 %) | 6.044 | 5.994 (0.827 %) |
PZT-5A | 0.0849 | 5.852 (3.305 %) | 6.052 | 5.994 (0.958 %) |
PZT-5H | 0.1751 | 5.741 (5.076 %) | 6.048 | 5.994 (0.893 %) |
G1195N | 0.2710 | 5.636 (6.812 %) | 6.048 | 5.994 (0.893 %) |
The convergence graphs plotted in Figures 3 and ^{4} for the tip deflection and potential developed at the root, respectively, compare the efficiency of the FSDT-based piezoelectric beam finite elements. The G1195N bimorph which has the highest modulus ratio among the chosen materials is taken as a particular example for this study. The FSDT-Coupled shows single-element convergence, closely reproducing the ANSYS-2D solutions for both the tip deflection and the potential developed. The conventional FSDT (^{Narayanan and Balamurugan, 2003}) overestimates the response as it neglects induced potential effects.
Modal Analysis: The accuracy and efficiency of the FSDT elements in predicting the natural frequencies of the bimorph cantilever beam shown in Figure 2 are compared here. The natural frequencies are evaluated for closed and open circuit electrical boundary conditions, with different materials. For open circuit, only the interface of the bimorph is grounded while, for closed circuit all the faces of bimorph are grounded. The results tabulated in Table 6 reveal the inability of conventional FSDT formulation to maintain the accuracy over the different bimorph materials. The consistent accuracy of the present FSDT-Coupled results validates the use of coupled polynomial shape functions in generating the element mass matrix consistent with the element stiffness matrix.
Materials | Modulus ratio | FSDT (Narayanan and Balamurugan, 2003) | ANSYS 2D | FSDT-Coupled |
---|---|---|---|---|
Open circuit | ||||
PVDF | 0.0100 | 169.74 (0.135 %) | 169.97 | 169.94 (0.018 %) |
PZT-2 | 0.0589 | 550.61 (0.736 %) | 554.69 | 554.41 (0.051 %) |
PZT-8 | 0.0687 | 556.12 (0.822 %) | 560.73 | 560.59 (0.025 %) |
PZT-4 | 0.0773 | 542.14 (0.852 %) | 546.80 | 547.38 (0.106 %) |
PZT-5A | 0.0849 | 463.68 (1.012 %) | 468.42 | 468.22 (0.043 %) |
PZT-5H | 0.1751 | 482.37 (1.814 %) | 491.28 | 491.50 (0.045 %) |
G1195N | 0.2710 | 506.16 (2.628 %) | 519.82 | 520.03 (0.040 %) |
Closed circuit | ||||
PVDF | 0.0100 | 169.11 (0.142 %) | 169.35 | 169.32 (0.018 %) |
PZT-2 | 0.0589 | 539.01 (0.739 %) | 543.02 | 542.90 (0.022 %) |
PZT-8 | 0.0687 | 542.48 (0.960 %) | 547.74 | 547.06 (0.124 %) |
PZT-4 | 0.0773 | 527.63 (0.991 %) | 532.91 | 532.63 (0.053 %) |
PZT-5A | 0.0849 | 449.75 (1.189 %) | 455.16 | 454.44 (0.158 %) |
PZT-5H | 0.1751 | 453.81 (2.270 %) | 464.35 | 463.53 (0.177 %) |
G1195N | 0.2710 | 461.90 (3.206 %) | 477.20 | 477.12 (0.017 %) |
Figures 5 and ^{6} show the comparison of convergence characteristics of FSDT-based piezoelectric beam finite element formulations to predict the first natural frequency of the G1195N bimorph in open and closed circuit conditions, respectively. FSDT-Coupled shows quick convergence, closely reproducing the ANSYS-2D solutions for both open and closed circuit conditions. The conventional FSDT (^{Narayanan and Balamurugan, 2003}) model underestimates the response as it neglects induced potential effects.
6.2 Example 2: A Two-Layer Asymmetric Piezoelectric Beam
A two-layer asymmetric piezoelectric cantilever beam having a steel host layer with a surface bonded piezoelectric layer of G1195N at the top, as shown in Figure 7 is considered here. The material properties used are:
Steel (^{Carrera and Brischetto, 2008}): E = 210GPa, υ = 0.3, ρ = 7850kgm ^{-3}
G1195N (^{Peng et al., 1998}): E = 63GPa, υ = 0.3, d _{31} = 254×10^{-12} mV ^{-1}, e _{31} = E×d _{31} Cm ^{-2}, ∈_{3} = 15×10^{-9} Fm ^{-1}, ρ = 7600kgm ^{-3}.
The length and total height of the beam are fixed (L = 100mm, h = 5mm), while thicknesses of the piezoelectric layer (h_{p} ) and the host layer (h_{c} ) are varied. The performance of the FSDT-based piezoelectric beam finite elements is evaluated over a wide range of the piezoelectric material proportion in the total beam thickness (thickness ratio:r = h_{p} / h). For a comparative evaluation of accuracy of various FSDT based formulations, converged results from a refined mesh of 40 elements have been used. The converged results from an ANSYS 2D simulation with a mesh of 200×20 elements are used as a benchmark.
Static Analysis-Sensor Configuration: In this configuration, the beam shown in Figure 7 is subjected to a tip load of -1000 N. The results for transverse tip deflection, axial tip deflection and potential developed across the piezoelectric layer at the mid-span of the beam for various thickness ratios are tabulated in Tables 7, ^{8} and ^{9}, respectively. The present FSDT-Coupled formulation proves its versatility, yielding consistently accurate predictions over the entire range of thickness ratio. The conventional FSDT formulation (^{Narayanan and Balamurugan, 2003}) does not maintain the consistent accuracy, as it neglects the induced potential coupling. The associated error increases rapidly in the higher thickness ratio regimes.
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
0.05 | -0.1671 | -0.1671 (0.000 %) | -0.1671 (0.000 %) |
0.1 | -0.1819 | -0.1819 (0.000 %) | -0.1819 (0.000 %) |
0.2 | -0.2109 | -0.2108 (0.047 %) | -0.2110 (0.047 %) |
0.3 | -0.2355 | -0.2354 (0.042 %) | -0.2362 (0.297 %) |
0.4 | -0.2523 | -0.2521 (0.079 %) | -0.2543 (0.793 %) |
0.5 | -0.2604 | -0.2601 (0.115 %) | -0.2647 (1.651 %) |
0.6 | -0.2622 | -0.2618 (0.153 %) | -0.2699 (2.937 %) |
0.7 | -0.2630 | -0.2625 (0.190 %) | -0.2757 (4.829 %) |
0.8 | -0.2704 | -0.2700 (0.148 %) | -0.2914 (7.766 %) |
0.9 | -0.2990 | -0.2987 (0.100 %) | -0.3378 (12.98 %) |
1.0 | -0.4008 | -0.4005 (0.075 %) | -0.5088 (26.95 %) |
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
0.05 | 0.1896 | 0.1897 (0.053 %) | 0.1897 (0.053 %) |
0.1 | 0.4040 | 0.4042 (0.050 %) | 0.4043 (0.074 %) |
0.2 | 0.8915 | 0.8918 (0.034 %) | 0.8926 (0.123 %) |
0.3 | 1.4065 | 1.4063 (0.014 %) | 1.4111 (0.327 %) |
0.4 | 1.8648 | 1.8631 (0.091 %) | 1.8792 (0.772 %) |
0.5 | 2.1898 | 2.1812 (0.393 %) | 2.2197 (1.365 %) |
0.6 | 2.3170 | 2.3148 (0.095 %) | 2.3867 (3.008 %) |
0.7 | 2.2523 | 2.2522 (0.004 %) | 2.3656 (5.030 %) |
0.8 | 1.9790 | 1.9803 (0.065 %) | 2.1376 (8.014 %) |
0.9 | 1.4033 | 1.4044 (0.078 %) | 1.5885 (13.20 %) |
1.0 | 0.0000 | 0.0000 (0.000 %) | 0.0000 (0.000 %) |
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
0.05 | 16.356 | 16.355 (0.006 %) | 16.355 (0.006 %) |
0.1 | 34.840 | 34.848 (0.023 %) | 34.851 (0.032 %) |
0.2 | 76.880 | 76.880 (0.000 %) | 76.950 (0.091 %) |
0.3 | 121.25 | 121.24 (0.008 %) | 121.66 (0.338 %) |
0.4 | 160.62 | 160.62 (0.000 %) | 162.01 (0.865 %) |
0.5 | 188.05 | 188.05 (0.000 %) | 191.36 (1.761 %) |
0.6 | 199.57 | 199.57 (0.000 %) | 205.76 (3.102 %) |
0.7 | 194.17 | 194.16 (0.005 %) | 203.94 (5.032 %) |
0.8 | 170.73 | 170.72 (0.006 %) | 184.29 (7.942 %) |
0.9 | 121.08 | 121.07 (0.008 %) | 136.95 (13.11 %) |
1.0 | 000.00 | 000.00 (0.000 %) | 000.00 (0.000 %) |
The convergence graphs plotted in Figures 8 and ^{9} for the transverse tip deflection and the potential developed at the root, respectively, prove the consistent efficiency of the present FSDT-Coupled formulation, which exhibits single-element convergence to ANSYS-2D solutions. FSDT (^{Narayanan and Balamurugan, 2003}) model shows very slow convergence to the inaccurate results, due to induced potential effects.
Static Analysis-Actuator Configuration: Here, the beam shown in Figure 7 is actuated by 100 volts. The variations of transverse and axial deflections at the tip, with thickness ratio are tabulated in Tables 10 and ^{11}, respectively. The FSDT-Coupled formulation consistently gives accurate predictions of results over the entire range of thickness ratio as given by ANSYS 2D simulation. The conventional FSDT formulation (^{Narayanan and Balamurugan, 2003}) does not yield consistently accurate results.
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
0.05 | 9.9834 | 9.9843 (0.009 %) | 9.9846 (0.012 %) |
0.1 | 10.839 | 10.838 (0.009 %) | 10.842 (0.028 %) |
0.2 | 12.463 | 12.463 (0.000 %) | 12.475 (0.096 %) |
0.3 | 13.673 | 13.674 (0.007 %) | 13.725 (0.380 %) |
0.4 | 14.099 | 14.101 (0.014 %) | 14.238 (0.986 %) |
0.5 | 13.557 | 13.558 (0.007 %) | 13.828 (1.999 %) |
0.6 | 12.156 | 12.156 (0.000 %) | 12.582 (3.504 %) |
0.7 | 10.180 | 10.180 (0.000 %) | 10.747 (5.570 %) |
0.8 | 7.8347 | 7.8349 (0.002 %) | 8.5019 (8.516 %) |
0.9 | 4.9635 | 4.9638 (0.006 %) | 5.6344 (13.52 %) |
1.0 | 0.0000 | 0.0000 (0.000 %) | 0.0000 (0.000 %) |
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
0.05 | -0.1752 | -0.1752 (0.000 %) | -0.1752 (0.000 %) |
0.1 | -0.2004 | -0.2005 (0.049 %) | -0.2005 (0.049 %) |
0.2 | -0.2583 | -0.2584 (0.039 %) | -0.2585 (0.077 %) |
0.3 | -0.3198 | -0.3201 (0.094 %) | -0.3206 (0.250 %) |
0.4 | -0.3763 | -0.3762 (0.027 %) | -0.3778 (0.399 %) |
0.5 | -0.4168 | -0.4170 (0.048 %) | -0.4206 (0.912 %) |
0.6 | -0.4378 | -0.4388 (0.228 %) | -0.4450 (1.645 %) |
0.7 | -0.4445 | -0.4455 (0.225 %) | -0.4537 (2.069 %) |
0.8 | -0.4459 | -0.4461 (0.045 %) | -0.4546 (1.951 %) |
0.9 | -0.4535 | -0.4542 (0.154 %) | -0.4599 (1.411 %) |
1.0 | -0.5080 | -0.5080 (0.000 %) | -0.5080 (0.000 %) |
Modal Analysis: The FSDT-Coupled formulation is evaluated here for its accuracy and efficiency to predict the natural frequencies of the asymmetric piezoelectric smart beam. The first natural frequency of the asymmetric Steel/G1195N beam shown in Figure 7 is computed for both open and closed circuit electrical boundary conditions. Table 12 shows the variation of first natural frequencies with thickness ratio. The results of FSDT-Coupled formulation agree very well with the ANSYS 2D simulation results. The results of the conventional FSDT formulation (^{Narayanan and Balamurugan, 2003}), show significant deviation in the higher thickness ratio regimes where the induced potential effects are predominant.
Thickness ratio (r) | ANSYS 2D | FSDT-Coupled | FSDT (Narayanan and Balamurugan, 2003) |
---|---|---|---|
Open circuit | |||
0.05 | 399.07 | 398.93 (0.035 %) | 398.93 (0.035 %) |
0.1 | 382.81 | 382.67 (0.037 %) | 382.66 (0.039 %) |
0.2 | 356.16 | 356.03 (0.037 %) | 355.87 (0.081 %) |
0.3 | 337.57 | 337.45 (0.036 %) | 336.89 (0.201 %) |
0.4 | 326.64 | 326.57 (0.021 %) | 325.17 (0.450 %) |
0.5 | 322.02 | 322.03 (0.003 %) | 319.24 (0.863 %) |
0.6 | 321.37 | 321.48 (0.034 %) | 316.61 (1.481 %) |
0.7 | 321.38 | 321.54 (0.050 %) | 313.76 (2.371 %) |
0.8 | 317.72 | 317.58 (0.044 %) | 305.71 (3.780 %) |
0.9 | 302.33 | 302.41 (0.026 %) | 284.40 (5.931 %) |
1.0 | 261.63 | 261.62 (0.003 %) | 232.14 (11.27 %) |
Closed circuit | |||
0.05 | 396.50 | 396.36 (0.035 %) | 396.36 (0.035 %) |
0.1 | 377.63 | 377.49 (0.037 %) | 377.47 (0.042 %) |
0.2 | 345.86 | 345.71 (0.043 %) | 345.55 (0.090 %) |
0.3 | 322.79 | 322.65 (0.043 %) | 322.06 (0.226 %) |
0.4 | 308.73 | 308.62 (0.035 %) | 307.13 (0.518 %) |
0.5 | 302.91 | 302.86 (0.016 %) | 299.88 (1.000 %) |
0.6 | 303.27 | 303.30 (0.010 %) | 298.13 (1.695 %) |
0.7 | 306.43 | 306.53 (0.033 %) | 298.35 (2.637 %) |
0.8 | 307.46 | 307.59 (0.042 %) | 295.31 (3.952 %) |
0.9 | 298.35 | 298.43 (0.027 %) | 280.17 (6.094 %) |
1.0 | 261.63 | 261.62 (0.003 %) | 232.14 (11.27 %) |
The consistent efficiency of the present FSDT-Coupled is revealed by the convergence graphs for first natural frequency in both open and closed circuit electrical boundary conditions plotted in Figures 10 and ^{11}, respectively. As seen from the figures, the FSDT-Coupled gives fast convergence, while FSDT (^{Narayanan and Balamurugan, 2003}) model shows very slow convergence to the inaccurate results, due to induced potential effects.
7 CONCLUSIONS
Based on coupled polynomial field interpolations in conjunction with a consistent through-thickness electric potential, a novel FSDT based extension mode piezoelectric beam finite element has been presented. The derived set of coupled shape functions handles bending-extension, bending-shear and induced potential couplings in a variationally consistent manner. Numerical evaluation has proven the merits of the present formulation over the conventional formulations available in the literature, in terms of accuracy and efficiency. From the numerical analysis, it was found that:
The performance of the conventional FSDT piezoelectric beam finite elements depends on the geometric and material parameters of the beam. The accuracy depends on the proportion of piezoelectric material in the total beam thickness (thickness ratio) and the ratio of induced modulus to the elastic modulus (modulus ratio). The convergence rate of conventional FSDT elements is deteriorated by the presence of bending-shear and bending-extension couplings.
The performance of the proposed FSDT-Coupled formulation proves to be insensitive to the material and geometric configuration of the beam cross-section. It consistently maintains the level of accuracy and efficiency for all modulus and thickness ratios.