Numerical analysis of end effects in laminated piezoelectric circular cylinders

Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186 www.elsevier.com/locate/cma

Numerical analysis of end effects in laminated piezoelectric circular cylinders C.W. Liu, E. Taciroglu

*

Department of Civil and Environmental Engineering, University of California, 5731E Boelter Hall, Box 951593, Los Angeles, CA 90095, USA Received 5 July 2005; received in revised form 27 October 2006; accepted 21 November 2006

Abstract A method to quantify Saint-Venant’s principle for laminated piezoelectric circular cylinders is described. An algebraic eigensystem is constructed based on homogeneous equations of equilibrium with displacement and electric potential (voltage) fields stated in exponential form. The end effects for displacements and voltages, as well as stresses and electric displacements of self-equilibrated states are represented by the eigendata extracted from the eigensystem. The real parts of the eigenvalues convey information on the inverse decay lengths and their corresponding eigenvectors are displacement and voltage distributions of the self-equilibrated states. Stress and electric displacement eigenvectors are formed by appropriate differentiation of the eigenvectors for the displacement and voltage fields and through the coupled electromechanical constitutive laws. The right and left-eigenvectors, which are obtained from the eigensystem and its adjoint, are related through bi-orthogonality relationships. The end effects due to arbitrary displacements, voltages, tractions or electric displacements prescribed at the cylinder’s ends are obtained by means of an expansion theorem based on these bi-orthogonality relationships, or by a least-squares method. A number of verification examples are provided to demonstrate that the present results compare well with known analytical solutions and numerical results obtained via three-dimensional finite element analyses.  2006 Elsevier B.V. All rights reserved. Keywords: Piezoelectricity; Saint-Venant’s principle; End effects

1. Introduction There have been extensive studies on the theoretical and numerical aspects of the Saint-Venant’s principle and related issues in elasticity theory. Pioneering works of Toupin [25] and Knowles [18] provided findings on the qualitative aspects of Saint-Venant’s principle through energy-based approaches. Johnson and Little [16] and Little and Childs [19] performed the first quantitative analyses of end effects in homogeneous, isotropic plane-strain strips and circular cylinders, and provided analytical solutions. Comprehensive reviews of this body of work may be found in Horgan and Knowles [15], and Horgan [13,14]. There

*

Corresponding author. Tel.: +1 310 2674655; fax: +1 310 2062222. E-mail address: [email protected] (E. Taciroglu).

0045-7825/$ - see front matter  2006 Elsevier B.V. All rights reserved. doi:10.1016/j.cma.2006.11.012

have also been numerical studies on Saint-Venant’s principle (see, for example, [10,12,17]). Studies on Saint-Venant’s principle for piezoelectricity, on the other hand, are comparatively scarcer. Notable examples include Batra and Yang [3] who extended Toupin’s approach to linear piezoelectricity; Fan [11] who investigated the decay rates in piezoelectric strips; and Ruan et al. [21] who examined the end effects in a twodimensional, semi-infinite piezoceramic strip by employing the Airy’s stress function approach. More recently, Tarn and Huang [24] developed a state-space approach for generalized plane-strain analyses of end effects and stress decay in laminated piezoelectric strips; and Borrelli et al. [5–7] studied anti-plane shear problems for linear, homogeneous piezoelectric solids. The quantification of Saint-Venant’s principle is especially pertinent for piezoelectric structures, as the stresses due to end effects penetrate into the interior further than those for purely elastic structures (see, for

2174

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

example, [9]). These stresses and their associated characteristic decay lengths are deemed important in the design and consideration of supports and geometry transition zones in piezoelectric sensors and actuators, and in the prediction of delamination failures. In Taciroglu et al. [22], a semi-analytic finite element method has been developed for analysis of multi-layered circular piezoelectric cylinders under general axisymmetric electromechanical loads. This method has significant computational and analytical advantages in the design of multilayered sensor/actuators and in passive/active control applications [23]. The said method has been constructed through an extension of the relaxed formulation of SaintVenant and Almansi-Michell problems, where the end conditions were given in terms of integral resultants of axial force, torque and electric flux. Under such restrictions, the solutions are rigorously valid everywhere except in the neighborhoods of the two ends. In order to obtain solutions for the whole cylinder, tractions and electric flux, prescribed pointwise at the ends, must be considered. The present work provides an analysis of these end effects; and as such, it brings the semi-analytic method for axisymmetric problems of Taciroglu et al. [22] to full-generality. The proposed method for determination of the end solutions is applicable for laminated piezoelectric cylinders where each layer may possess the most general form of anisotropy, unlike any prior analytical treatments mentioned above whereby the solutions are limited to special cases of material anisotropy and cross-sectional geometry. Additionally, the semi-analytical method is computationally more efficient than fully discrete methods, because part of the solution (i.e., along the axis of the cylinder) is obtained analytically, and the problem geometry is discretized only along the radial direction. Finally, the proposed method yields direct measures of decay rates of the end effects, unlike fully discrete approaches which require post-processing of the determined stress fields. In what follows, an exponential form of the field variables is substituted into the homogeneous form of governing equations in Taciroglu et al. [22] to obtain an algebraic eigensystem for a quantitative analysis of end effects. As such, an examination of Taciroglu et al. [22] should help the reader to establish the prerequisites for the work presented here. The solutions of the eigensystem contain the inverse decay lengths and their corresponding eigenvectors. This eigendata is then used in an expansion theorem to construct the end solutions. Results of several verification studies are provided to illustrate the use and the accuracy of the proposed approach. 2. Preliminaries Consider a laminated piezoelectric cylinder of length L as shown in Fig. 1. Self-equilibrated tractions and electric flux at the tip-end (z = 0), and displacements and voltage at the root-end (z = L) are prescribed pointwise. The cylinder may be composed of any number of perfectly bonded

Fig. 1. Problem geometry and the coordinate system.

layers, each possessing a distinct thickness and linear piezoelectric properties. Cylindrical coordinates (r, h, z) are used with the origin located at the center of the cross-section on the tip-end, and the z-axis running toward the root-end. The primary variables are the mechanical displacements T u ¼ ½u; v; w , and the electric potential /. The dependent T variables are the electric displacement D ¼ ½Dr ; Dh ; Dz  , T electric field E ¼ ½Er ; Eh ; Ez  ¼ r/, stress r ¼ ½rrr ; rhh ; rzz ; rhz ; rrz ; rrh T , and the strain e ¼ ½err ; ehh ; ezz ; ehz ; erz ; erh T . The constitutive equations of linear piezoelectricity are stated in the following general form: r ¼ ce  eE; D ¼ eT e þ jE;

ð1Þ

where c, e and j are the matrices of the elastic anisotropic moduli (6 · 6), piezoelectric constants (3 · 6) and dielectric constants (3 · 3), respectively. By combining the mechanical and electrical variables as       e r u q¼ ; Q¼ ; v¼ ; ð2Þ E 91 D 91 / 41 the constitutive relation for a given lamina of the piezoelectric cylinder can be rewritten as   c e   Q ¼ C q where C ¼ T : ð3Þ e j In addition, there are nine generalized deformational relations, q ¼ Lfvg, where the linear operator L performs the appropriate differentiations in cylindrical coordinates that relate the strain and electric displacement fields to the mechanical displacement and potential fields. For piezoelectric materials, large discrepancies between the numerical values of c, e and j exists, and this typically causes numerical anomalies in finite-precision arithmetic (see, for example,[2]). The use of dimensionless (normalized) variables precludes such effects as ill-conditioning encountered in least-squares solutions which are utilized in the present work. In setting forth this normalization, three key properties are selected as the reference values, viz., (1) total cylinder thickness t, (2) an elastic modulus c0, and (3) a piezoelectric constant e0, where c0 and e0 are of a particular lamina in the cylinder’s radial profile. The geometry, mechanical displacements, electric potential, and the material constants are normalized as

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

r z ui ¼ / ; r ¼ ; z ¼ ;  ui ¼ ði ¼ r;h;zÞ; / ð4Þ t t t E0 t cpq jij eip ij ¼ 0 ; eip ¼ 0 ðp;q ¼ 1;2; ...; 6; i; j ¼ 1;2; 3Þ; cpq ¼ 0 ; j c j e ð5Þ where E0 ¼ c0 =e0 , and j0 is the reference dielectric constant 2 given by j0 ¼ ðe0 Þ =c0 . The remaining variables are rendered dimensionless by rp Dk ; ep ¼ ep ; Dk ¼ 0 ; c0 e Ek Ek ¼ 0 ðp ¼ 1; 2; . . . ; 6; k ¼ 1; 2; 3Þ: E p ¼ r

ð6Þ

The overbar, which identifies normalized variables, is dropped for the sake of brevity. This scheme yields completely dimensionless equations in the same form as their dimensional counterparts. 3. Algebraic eigensystem The eigendata for representing the end effects are determined through the homogeneous form of the governing kinematic equations of equilibrium, which were derived in Taciroglu et al. [22] as K1 V þ K3 V;z  K6 V;zz ¼ 0;

ð7Þ

where V(z) is the vector of nodal displacements and voltages. These equations are based on a semi-analytical finite element formulation, where the discretization of the laminated cylinder takes the form of a series of three-node cylindrical laminae with quadratic interpolation functions. The Ki’s are stiffness matrices, where K1 and K6 are symmetric, while K3 is antisymmetric. The solution to Eq. (7) for a self-equilibrated state (i.e., a state that satisfies the homogeneous equations) has the form VðzÞ ¼ V0 ecz ;

ð8Þ

where V0 is the array representing the nodal displacements and electric potential distribution in the radial direction and c is the inverse decay length. Substituting the expression for V(z) in Eq. (8) into Eq. (7), and eliminating the common factor ecz yields a second-order algebraic eigenvalue problem 2

ðK1  cK3  c K6 ÞV0 ¼ 0:

ð9Þ

This quadratic form can be linearized and reduced to a first-order algebraic eigenvalue problem by introducing a new vector V1, which satisfies V1 ¼ cV0 ) K1 V1 ¼ cK1 V0 :

ð10Þ

Plugging Eq. (10) into Eq. (9) yields K1 V0  K3 V1  cK6 V1 ¼ 0:

ð11Þ

Combining Eqs. (10) and (11) yields a first-order algebraic eigenvalue problem



0

K1

K1 K3



V0 V1



 ¼c

K1

0

0

K6



2175

V0 V1

 or KUr ¼ cMUr ; ð12Þ

where Ur is the right-handed generalized coordinate vector. For a discretization with N degrees-of-freedom, 2N eigenvalues are contained in the algebraic eigensystem. These can be real, complex conjugate pairs or zero-valued. Zero roots represent constant strain modes and rigid-body motions, i.e., axial translation and rotation about the z-axis as well as constant potential throughout the cylinder. Nonzero roots represent attenuation rates of self-equilibrated end effects into the interior of the cylinder. Specifically, the real and complex roots respectively correspond to monotonic and oscillatory decays. Moreover, positive real roots and complex roots with positive real parts correspond to decay from the tip-end (z = 0), where a self-equilibrated stress and electric displacement state is applied. The negative real roots and the complex roots with negative real parts correspond to decay from the root-end (z = L), where a self-equilibrated displacement and electric potential state is applied. The eigenvalue with the smallest magnitude real part defines the inverse decay length with the furthest penetration into the interior. We rearrange the eigendata such that the real parts of eigenvalues increase from the smallest to the largest. Accordingly, their corresponding eigenvectors (U) are also sorted. The eigenvectors may be separated into upper and lower parts, denoted by Uu and Ul, which correspond to V0 and V1, respectively. Furthermore, each part may be divided into three portions as, for example, Uu ¼ ½Uuroot UuSV Uutip . Uuroot and Uutip contain the eigenvectors for displacement and electric potential, which will be used in the analysis of end effects, and each has the dimension of N · (N  3). The vector Uuroot corresponds to the negative real eigenvalues and complex eigenvalues with negative real parts. Similarly, the vector Uutip corresponds to the positive real eigenvalues and complex eigenvalues with positive real parts. The N · 6 matrix UuSV , which corresponds to the zero eigenvalues, contains the rigid body and constant strain modes. Lastly, the following relationships between the upper and lower eigenvectors can be obtained from Eq. (10) Ulroot ¼ Uuroot C1 ;

Ultip ¼ Uutip C2 ;

ð13Þ

where C1 and C2 are diagonal matrices of negative and positive eigenvalues, respectively. Consider the adjoint problem to Eq. (9) UT0 ðK1  cK3  c2 K6 Þ ¼ 0:

ð14Þ

Same linearization can be used here, so that Eq. (14) is reduced to a first-order algebraic eigenvalue problem       K1 0 U0 0 K1 U0 ¼c or KT Ul ¼ cMUl ; 0 K6 U1 K1 K3 U1 ð15Þ

2176

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

where Ul is the left-handed generalized coordinate vector. The eigenvalues of Eq. (15) are the same as those of Eq. (12) with W as the left modal (eigenvector) matrix. Same partitions may be made as those on U and thus, Wlroot ¼ Wuroot C1 , Wltip ¼ Wutip C2 hold. The right and left eigenvectors satisfy the standard biorthogonality relations WT MU ¼ I

ð16Þ

from which the bi-orthogonality conditions expressed in terms of the upper half eigendata may be stated as Wuroot T K3 Uuroot þ cm cn Wuroot T K1 Uuroot ¼ dmn ; Wutip T K3 Uutip þ cm cn Wutip T K1 Uutip ¼ dmn :

ð17Þ

A dual expansion theorem can be formulated with these bi-orthogonality relations. To wit, a given arbitrary displacement vector f can be expressed in terms of the righthanded eigenvectors Uuroot in the form f ¼ Uuroot a;

ð18Þ

where a represents the amplitudes. These amplitudes may be obtained directly through the bi-orthogonality conditions as a ¼ ðWuroot T K3 þ c2n Wuroot T K1 Þf:

ð19Þ

4. Analysis of end effects V(z) denotes the mechanical displacement and electric potential fields, and Q(z) denotes the stress and electrical displacement fields. As mentioned previously, V(z) contains all the nodal mechanical displacements and electric potentials of the discretized model. Similarly, Q(z) is a combined array of (six) stress and (three) electrical displacement components at each Gaussian quadrature point, taken over all such points of an element and over all elements comprising the discretized model. In the present approach (and without any loss of generality), V(z) and Q(z) are subject to boundary conditions prescribed pointwise at the tip and the root-ends, respectively. Due to linearity, the solution to a given problem with arbitrary pointwise specification of the end conditions can be stated as the superposition of Saint-Venant and end solutions, VðzÞ ¼ VSV ðzÞ þ Vend ðzÞ;

QðzÞ ¼ QSV ðzÞ þ Qend ðzÞ:

ð20Þ

Here, the subscript SV denotes the Saint-Venant solution, which may be obtained via the semi-analytical method described in Taciroglu et al. [22]; whereas the subscript end represents the end solution due to prescribed selfequilibrated tractions and restraints, which is the subject matter of the present work. It is noted that the term ‘‘end solution’’ denotes the difference between the Saint-Venant solution and that for the arbitrary and not necessarily self-equilibrated end conditions prescribed pointwise. As

such, the end solutions contain the response due to purely self-equilibrated boundary conditions at the ends of the cylinder. The end solutions may be separated into two parts, one for decay from the tip-end and the other for decay from the root-end, which may be expressed as     atip cj z cj ðLzÞ ; VðzÞ ¼ VSV ðzÞ þ Utip diagðe Þ Uroot diagðe Þ aroot ð21Þ     atip cj z cj ðLzÞ ; QðzÞ ¼ QSV ðzÞ þ Rtip diagðe Þ Rroot diagðe Þ aroot ð22Þ where Uroot and Utip are identical with Uuroot and Uutip , respectively, and atip and aroot are the unknown amplitudes. Rtip and Rroot contain the modal stress and electrical displacement of the self-equilibrated states, which are obtained from the modal displacement and electric potential at the element level, Rtip ¼ C ðb1 Utip þ b2 Utip ;z Þ ¼ C ðb1 Utip  cb2 Utip Þ; Rroot ¼ C ðb1 Uroot þ b2 Uroot;z Þ ¼ C ðb1 Uroot  cb2 Uroot Þ: ð23Þ Details of this formulation and the transformation matrices b1 and b2 can be found in Taciroglu et al. [22]. According to Saint-Venant’s principle, the end effects are reasonably confined to the neighborhood of each end. Here, we first assume that the cylinder sufficiently long, such that it is possible to consider the effects on one end completely uncoupled from the other end. 4.1. Displacement and electric potential boundary conditions at the root-end (z = L) For a sufficiently long cylinder, Eq. (21) can be restated as VðzÞ ¼ VSV ðzÞ þ Uroot diagðecj ðLzÞ Þaroot :

ð24Þ

Let VL denote the nodal values of the arbitrary displacement and electric potential restrains, which are prescribed at the root-end, we can obtain the root-end condition by evaluating Eq. (24) at z = L. To wit, VL ¼ VSV ðLÞ þ Uroot aroot or Uroot aroot ¼ VL  VSV ðLÞ:

ð25Þ

Unlike the stress and electric displacement conditions, VSV(L) cannot be obtained easily. However, it can be expressed in the form VSV ðLÞ ¼ USV b;

ð26Þ

where USV is the rigid body and constant strain modes and b are the unknown coefficients. USV contains the eigenvectors corresponding to zero eigenvalues of the quadratic algebraic eigenvalue problem stated in Eq. (9), from which the following condition can be obtained

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

K1 USV ¼ 0:

ð27Þ

Pre-multiplying both sides of Eq. (25) with K1 yields K1 Uroot aroot ¼ K1 VL :

ð28Þ

A least-squares method can be applied to Eq. (28) in order to obtain the amplitudes aroot . To wit, pre-multiplying both sides of Eq. (28) with the conjugate transpose of H K1 Uroot , denoted by ½K1 Uroot  and then inverting the resulting matrix yields 1

aroot ¼ ½UHroot KH1 K1 Uroot  UHroot KH1 K1 VL :

ð29Þ

If only self-equilibrated displacement and electric potential restrains on VL are prescribed at the root-end, then the root-end condition Eq. (25) becomes VL ¼ Uroot aroot :

ð30Þ

In this case, the amplitudes aroot can be obtained directly by making use of the bi-orthogonality conditions in Eq. (17) as aroot ¼ ðWHroot K3 þ c2n WHroot K1 ÞVL :

ð31Þ

Alternatively, the amplitudes may also be obtained by a least-squares method as aroot ¼ ½UHroot Uroot 1 UHroot VL :

ð32Þ

4.3. Mixed boundary conditions Here we consider arbitrary mixed boundary conditions prescribed at the tip-end (z = 0). We assume that displacement and electric potential conditions, denoted by V10 , are enforced at a set of finite element nodes; stress and electric displacement conditions, denoted by S20 , are prescribed at the Gaussian quadrature points within the rest of the finite elements. The uniform displacement and electric potential state VSV ð0Þ, the uniform stress and electric displacement state SSV ð0Þ are extracted from V10 and S20 , respectively, such that their resultants are the same as those of V10 and S20 . Then, the relative modal displacement and electric potential at nodes, denoted by U1tip , are composed by extracting the corresponding eigenvalues from Utip . The relative modal stress and electrical displacement at the Gaussian points, denoted by S2tip , are composed by extracting the corresponding eigenvalues from Stip . Subsequently, the mixed boundary conditions at the tip-end can be stated as the following: " #   " 1 # Utip VSV ð0Þ V10 ð38Þ ¼ atip ; þ 2 SSV ð0Þ S0 S2tip and thus, "

4.2. Traction and electric displacement conditions at the tip-end (z = 0) The stress and electric displacement representations can be obtained from Eq. (22), which has the form QðzÞ ¼ QSV ðzÞ þ Rtip diagðecj z Þatip :

ð33Þ

Only the stress components ðrrz ; rhz ; rzz Þ and the electric displacement component Dz may be prescribed at the tipend (z = 0). Let S0 denote these prescribed components’ values at all of the N Gaussian points of the discretized cross section, we have ST0 ¼ ½rrz1 ; rhz1 ; rzz1 ; Dz1 ; . . . rrzN ; rhzN ; rzzN ; DzN :

ð34Þ

In order to obtain the tip-end condition, the same components of QSV and Rtip need to be extracted, which are denoted by SSV and Stip . Subsequently, Eq. (33) becomes SðzÞ ¼ SSV ðzÞ þ Stip diagðecj z Þatip :

ð35Þ

Then, Eq. (35) is evaluated at z = 0, which leads to S0 ¼ SSV ð0Þ þ Stip atip

or

Stip atip ¼ S0  SSV ð0Þ:

ð36Þ

It should be noted that it is quite straightforward to separate the arbitrary S0 into a uniform stress and electric displacement state, SSV ð0Þ, with the same resultants as those of S0 and a self-equilibrated stress and electric displacement state with zero resultants. Performing a least-squares solution to Eq. (36) leads to 1

atip ¼ ½SHtip Stip  SHtip ½S0  SSV ð0Þ:

ð37Þ

2177

U1tip S2tip

#

" atip ¼

V10 S20

#

 

VSV ð0Þ SSV ð0Þ

 ð39Þ

:

As before, we can obtain the amplitudes atip by applying a least-squares solution to Eq. (39). 4.4. Coupled end effects Here, we consider a cylinder for which stress and electric displacement conditions S0, and displacement and electric potential conditions VL are applied at the tip and rootends, respectively; and that the cylinder is sufficiently short so that the two end effects are coupled. In this case, the coefficients atip and aroot must be solved simultaneously using Eqs. (21) and (22), which should be evaluated at z = L and z = 0, respectively. To wit, 

VL S0









Utip diagðecj L Þ Uroot þ ¼ SSV ð0Þ Stip Sroot diagðecj L Þ VSV ðLÞ



atip



aroot ð40Þ

and thus,        VL VSV ðLÞ atip Utip diagðecj L Þ Uroot ¼  : Stip Sroot diagðecj L Þ aroot S0 SSV ð0Þ ð41Þ The least-squares solution of Eq. (41) yields the coefficients atip and aroot . Consequently, Eqs. (21) and (22) yield the solutions of the displacement and electric potential as well as stress and electric displacement states for the whole cylinder.

2178

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

It is noted that U, R, and S are rectangular matrices whose column dimension is equal to the number of eigenvectors used. The accuracy of the end solutions depends on this number. Usually, a small number of lower modes, rather than all, can be used in the analysis without affecting the solution accuracy. However, it is observed that the tractions and electric displacements prescribed at the Gaussian quadrature points demand a greater number of modes than the displacement and electric potential boundary conditions prescribed at the nodes. 5. Verification studies In this section, a purely mechanical verification example for a solid isotropic cylinder is considered (Section 5.1) first, and the results are compared with known analytical solutions. Then, decay rates of electric and mechanical fields for a transversely isotropic piezoelectric cylinder are obtained and compared with analytically obtained results (Section 5.2). Finally, the end effects for a homogeneous (Section 5.3), a laminated hollow cylinder (Section 5.4– 5.6), both made of lead–zirconium–titanate, and a hybrid laminated cylinder made of soft elastic material and lead– zirconium–titanate (Section 5.7) are analyzed with the present method. These results are compared with threedimensional numerical solutions, obtained using the proprietary finite element analysis package ANSYS [1]. 5.1. Stress boundary problem for an isotropic elastic solid cylinder Consider a homogeneous, solid, circular cylinder with radius ro ¼ 1. Without losing generality, we choose the Young’s modulus as E = 1.0 GPa and the Poisson’s ratio as m = 0.0 (Case 1), m = 0.3 (Case 2). First, we calculate the eigenvalues for both cases, which define the inverse decay lengths from the tip-end. A subset of the lowest eigenvalues for these cases are tabulated and along with analytical solutions by Little and Childs [19] in Table 1, and a very good agreement is observed between the two. In this analysis, 15 finite elements were used through the thickness. It should be noted here that the modes found in the present work (i.e., modes 3, 6, 9, etc. in Table 1), but are missing in Little and Childs [19], are the circumferential modes. These modes represent constant circumferential displacement states, which are absent in a torsion/ twist-free axisymmetric problem such as that studied by Little and Childs [19]. Subsequently, we consider a stress boundary problem where the components rzz ¼ 1  2r2 and rrz ¼ 2:4r 2:6r3 þ 0:2r5 are prescribed at the tip-end, which is the same problem examined in Little and Childs [19]. These particular radial variations were chosen by Little and Childs, because they are fairly general in form, and both components conveniently represent a purely self-equilibrated stress state. Figs. 2 and 3 display the variation of the normalized non-zero stress and displacement compo-

Table 1 Subsets of eigenvalues for decay from the tip-end Mode

Present method

Little and Childs

Case 1 1,2 3 4,5 6 7,8 9 10,11 12 13,14

2.55677 ± i1.38897 5.13565 6.00596 ± i1.63871 8.41762 9.23402 ± i1.82896 11.62192 12.42197 ± i1.96651 14.80323 15.59967 ± i2.06929

2.55677 ± i1.38897 – 6.00586 ± i1.63870 – 9.23317 ± i1.82906 – 12.41789 ± i1.96788 – 15.58596 ± i2.07680

Case 2 1,2 3 4,5 6 7,8 9 10,11 12 13,14

2.72218 ± i1.36220 5.13565 6.06025 ± i1.63751 8.41762 9.26829 ± i1.82676 11.62192 12.44980 ± i1.95770 14.80323 15.63121 ± i2.03573

2.72218 ± i1.36210 – 6.06008 ± i1.63762 – 9.26684 ± i1.82826 – 12.44253 ± i1.96724 – 15.60544 ± i2.07628

nents, respectively, through the thickness of the cylinder at cross-sections with different axial coordinates as computed via the present method and the analytical method of Little and Childs [19]. The agreement between the results of the two methods are remarkably good. The results in Figs. 2 and 3 are normalized with respect to the absolute maximum values of each stress and displacement component. These maximum values are omitted for brevity as the this example is only meant to display the veracity of the present approach. 5.2. The characteristic decay rates for transversely isotropic piezoelectric cylinders under circumferential shear The elastostatic problem considered in Section 5.1 involves only isotropic solids. An analytic solution for characteristic decay rates of end effects for transversely isotropic, semi-infinite piezoelectric strips subjected to antiplane shear has recently been obtained by Borrelli et al. [7] (see, also [8] for decay rates in functionally graded piezoelectric strips). For such problems, Borrelli et al. [7] have shown that the mechanical and electrical fields are uncoupled and that both fields have the same characteristic decay rates of p. Specifically, the mechanical decay rate cm =t and the electrical decay rate ce =t are equal to p (in the absence of pre-stress), where cm is the smallest eigenvalue corresponding to circumferential mode; ce is the smallest eigenvalue corresponding to electrical mode, and t is the thickness of the semi-infinite strip. Here, in order to validate the proposed semi-analytical method, we reproduce these results on transversely isotropic piezoelectric hallow circular cylinders, with the premise that as the inner radius to thickness ratio ri =t ! 1, the cylinder under a self-equilibrated circumferential shear load

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

2179

0

0.5

σ

σ

rr

θθ

0 –0.5

–0.5

–1

0

0.5 r

–1 0

1

1

0.5 r

1

0.5 r

1

1

0

0.5

σrz

σ

zz

0.5

–0.5 –1

0

0.5 r

0

1

0

Fig. 2. Verification for stress boundary problem of a solid cylinder. The results are the distributions of non-zero stress components ðrrr ; rhh ; rzz and rrz Þ over the thickness at different cross-sections (z/ri = 0 (s); 0.3 (}); 0.6 (,)). The symbols and the (dashed, dotted, etc.) lines denote the present solutions and the analytical solutions provided by Little and Childs [19], respectively.

0

0.5

w

u

0 –0.5

–0.5

–1 0

0.5

–1 0

1

0.5

r

1

r

Fig. 3. Verification for stress boundary problem of a solid cylinder. The results are the distributions of non-zero displacement components (u and w) over the thickness at different cross-sections (z/ri = 0 (s); 0.3 (}); 0.6 (,)). The symbols and the (dashed, dotted, etc.) lines denote the present solutions and the analytical solutions provided by Little and Childs [19], respectively.

approximates the behavior of a semi-infinite strip subject to anti-plane shear. For this verification problem, we choose the cylinder to be made of PZT-5H, a transversely isotropic piezoelectric material, with inner radius ri and thickness t. Table 2 displays the original and normalized material properties of PZT-5H, using the key reference parameters c0 ¼ c44 ¼ 23:0 GPa, e0 ¼ e33 ¼ 23:3 C=m2 , j0 ¼ 2:3604 108 F=m, and E0 ¼ 9:8712  108 N=C. Table 3 contains the mechanical and electrical decay rates obtained via the semi-analytical method for different inner radius to thickness ratios. These decay rates indeed approach to p as the said aspect ratio becomes very large. The results (i.e., decay rates normalized with p) displayed in Table 3 with asterisk and circle (*,s) marks are obtained with 30 finite elements used in the semi-analytical model, while the unmarked ones are obtained with 15 elements. Additionally, a purely mechanical problem (i.e., a trans-

Table 2 Original and normalized material constants for PZT-5H c11 c33 c44 c12 c13 c66 c22 c55

¼ 12:60  1010 (Pa) ¼ 11:70  1010 ¼ 2:30  1010 ¼ 7:95  1010 ¼ 8:41  1010 ¼ 2:33  1010 ¼ c11 ; c23 ¼ c13 ¼ c44

e13 e33 e15 e23 e42

¼ 6:50 (C/m2) ¼ 23:30 ¼ 17:00 ¼ e13 ¼ e15

j11 ¼ 15:03  109 (F/m) j33 ¼ 13:00  109 j22 ¼ j11

c11 c33 c44 c12 c13 c66 c22 c55

¼ 5:4783 ¼ 5:0870 ¼ 1:0 ¼ 3:4565 ¼ 3:6565 ¼ 1:0130 ¼ c11 ; c23 ¼ c13 ¼ c44

e13 e33 e15 e23 e42

¼ 0:2790 ¼ 1:0 ¼ 0:7296 ¼ e13 ¼ e15

j11 ¼ 0:6368 j33 ¼ 0:5508 j22 ¼ j11

2180

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

Table 3 Normalized mechanical and electrical decay rates of PZT-5H cylinders logðri =tÞ

Mechanical decay rate/p

Electrical decay rate/p

0 1 2 3 3 3

1.0517 1.0180 1.0052 1.0049 1.0049(*) 1.0049()

1.0899 1.0355 1.0329 1.0329 1.0329(*) –

versely isotropic cylinder with the same mechanical constitutive parameters as PZT-5H) is also analyzed with the present method and the mechanical decay rate is obtained (marked with a circle in Table 3). It can be seen that the mechanical decay rate under circumferential shear load is uncoupled from the electrical (potential) field. As such, all of the results obtained with the semi-analytical method are consistent with the analytical solutions of the class of problems considered in Borrelli et al. [7]. 5.3. Stress boundary problem for a hollow PZT cylinder Here, we consider a homogeneous hollow circular cylinder with length L = 20, and with inner and outer radii ri = 1, ro = 2, which are rendered dimensionless by a reference parameter t = 1 cm. The cylinder is made of PZT4 (lead–zirconium–titanate) material, which is a widely used piezoceramic in industrial applications. Table 4 displays the original [4] and normalized material properties of PZT4, using the key reference parameters c0 ¼ c44 ¼ 25:6 GPa, e0 ¼ e33 ¼ 15:1 C=m2 , j0 ¼ 8:90664  109 F=m, and E0 ¼ 1:69536  109 N=C. A self-equilibrated traction rzz ¼ 9r þ 14, ð1 6 r 6 2Þ, rendered dimensionless with c0 ¼ 25:6 GPa, is applied at the tip-end of the cylinder. To the best of the authors’ knowledge, analytical solution of this problem does not exist. So the present results are compared with those obtained via the three-dimensional finite element method [1]. In order to verify the present method with the normalization technique, original properties of PZT4 are used in Table 4 Original and normalized material constants for PZT4 c11 c33 c44 c12 c13 c66 c22 c55

¼ 13:90  1010 (Pa) ¼ 11:54  1010 ¼ 2:56  1010 ¼ 7:78  1010 ¼ 7:43  1010 ¼ 3:06  1010 ¼ c11 ; c23 ¼ c13 ¼ c44

e13 e33 e15 e23 e42

¼ 5:20 (C/m2) ¼ 15:10 ¼ 12:70 ¼ e13 ¼ e15

j11 ¼ 6:46  109 (F/m) j33 ¼ 5:62  109 j22 ¼ j11

c11 c33 c44 c12 c13 c66 c22 c55

¼ 5:42969 ¼ 4:49219 ¼ 1:0 ¼ 3:03906 ¼ 2:90234 ¼ 1:19531 ¼ c11 ; c23 ¼ c13 ¼ c44

e13 e33 e15 e23 e42

¼ 0:34437 ¼ 1:0 ¼ 0:84106 ¼ e13 ¼ e15

j11 ¼ 0:72530 j33 ¼ 0:63099 j22 ¼ j11

ANSYS model, where traction rzz ¼ c0 ð9r þ 14Þ; 1 6 r 6 2 is applied at the tip-end. As a consequence, the two sets of solutions have the following relationships: (i) ð1Þ ð2Þ uk ¼ uk ; ðk ¼ r; h; zÞ; (ii) E0 /ð1Þ ¼ /ð2Þ ; and (iii) c0 rpð1Þ ¼ ð2Þ rp ðp ¼ 1; 2; . . . ; 6Þ, where the superscripts (1) and (2) denote the present and the ANSYS solutions, respectively. The variation of the normalized non-zero stresses, displacements, and electric potential through the thickness of the cylinder at different cross-sections are displayed in Figs. 4 and 5, as computed via the two methods. The agreement between these two results are quite good. It is also noted that the present solutions yield rrr ¼ 0 at the inner and outer surfaces of the cylinder, and rrz ¼ 0 at its tip-end. This is in contrast with the ANSYS results. It is possible that the behavior predicted by ANSYS near the traction-free lateral surface may be more accurate, and indicate the presence of a boundary layer. Similar observations were made by Vel and Batra [26,27]. This issue will be left out of the scope of the present work; and further studies are needed to conclude regarding the presence/absence of boundary layers. 5.4. Decomposition procedure for arbitrary displacement and electric potential boundary conditions Prior to analyzing the displacement and electric potential boundary problem, we illustrate the claim, mentioned in Section 4 for decomposing arbitrary displacement and electric potential states into uniform and self-equilibrated states with zero resultants. It is noted that the solution for the uniform states (i.e., Saint-Venant solutions) can be obtained with methods described in Taciroglu et al. [22], whereas the present work addresses the solutions for self-equilibrated states, bringing the semi-analytical method to full generality. In order to illustrate the validity of the decomposition procedure, we first analyze the same hollow cylinder of the previous section (Section 5.3), subjected to normalized transverse displacements wL ðrÞ ¼ r2 þ 3r  2 applied at the cylinder’s root-end. We can obtain the corresponding self-equilibrated state from Eq. (29), and subsequently determine the uniform state that has the same resultant as the prescribed displacement state, wL(r). The present method’s results and those obtained via ANSYS are displayed in Fig. 6, and a comparison of the two validates the decomposition approach. Second, we consider a normalized electric potential /L ðrÞ ¼ r  1 applied at the root-end, and separate it into uniform and self-equilibrated parts via the proposed method. Fig. 7 displays the results of the present approach as well as those obtained via ANSYS. Again, the two sets of results agree with each other quite well. 5.5. Displacement boundary problem for a laminated PZT cylinder A piezoelectric hollow cylinder of two equal-thickness laminae with inner and outer radii ri ¼ 1, ro ¼ 2 is consid-

1.5

1.5

1

1 σθ θ

σrr

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

0.5

2181

0.5

0

0

–0.5 1

1.5 r

–0.5 1

2

1.5

1.5 r

2

1.5 r

2

1

1 0.5 σrz

σzz

0.5 0

0

–0.5 –1 1

1.5 r

–0.5 1

2

Fig. 4. Verification for stress boundary problem of a hollow PZT4 cylinder. The results are the distributions of non-zero stress components ðrrr ; rhh ; rzz and rrz Þ over the thickness at different cross-sections (z/ri = 0 (); 0.5 (); 1.0 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.

1

1

1

0.5

0.5

0

0

φ

w

u

0.5

0 –0.5 –0.5 1

1.5 r

–1 1

2

–0.5

1.5 r

–1 1

2

1.5 r

2

Fig. 5. Verification for stress boundary problem of a hollow PZT4 cylinder. The results are the distributions of non-zero displacement components (u and w) and electric potential (/) over the thickness at different cross-sections (z/ri = 0 (); 0.5 (); 1.0 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.

2

=

1.5

1 0

0.2

0.4

2

2

1.5

+ 1.5

1

0

0.1

0.2

1 –0.2

0

0.2

Fig. 6. Verification for displacement separation of a hollow PZT4 cylinder: (a) applied displacement state wL; (b) uniform displacement state; (c) selfequilibrated displacement state. In (b) and (c), the lines and the symbols denote the present and ANSYS results, respectively. The x and y-axes denote the displacement magnitude (w) and the radial distance (r), respectively.

ered under a self-equilibrated displacement state applied at the root-end. The inside and outside layers are made of

PZT4 material with crystal orientation angles hi ¼ 00 and ho ¼ 900 , respectively. Without any loss of generality, the

2182

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

2

2

2

1.5

= 1.5

+ 1.5

1 0

0.5

1

1 0

0.5

1

1 –1

–0.5

0

0.5

Fig. 7. Verification for electric potential separation of a hollow PZT4 cylinder: (a) applied electric potential state /L; (b) uniform electric potential state; (c) self-equilibrated electric potential state. In (b) and (c), the lines and the symbols denote the present and ANSYS results, respectively. The x and y-axes denote the electric potential magnitude (/) and the radial distance (r), respectively.

self-equilibrated displacement state is chosen to be one of the calculated eigenvectors. Again, the present results are compared with those obtained with ANSYS because an analytical solution of this problem does not exist. The variation of the normalized non-zero stresses, displacements and electric potential through the thickness of the cylinder at cross-sections with different axial coordinates are displayed in Figs. 8 and 9, as computed via the two methods. Again, a remarkably good agreement between the two sets of results is observed. 5.6. Electric potential boundary problem for a laminated PZT cylinder

age state, /self ðrÞ ¼ r2 þ 3r  2:166117, is extracted from a more general voltage state, /ðrÞ ¼ r2 þ 3r  2, via the decomposition technique described earlier. The present results are compared with those obtained with ANSYS, because an analytical solution of this problem does not exist. The variation of the normalized non-zero stresses, displacements, and the electric potential through the thickness of the cylinder at cross-sections with different axial coordinates are displayed in Figs. 10 and 11, as computed via the two methods. Again a good agreement between the two sets of results is observed. 5.7. Decay length investigation for a hybrid laminated cylinder

The same two-layer cylinder of the previous example is considered under a self-equilibrated electric potential state applied at the root-end. This applied self-equilibrated volt-

Consider a two-layered cylinder with the same geometry as the previous example. The inner and outer layers are

0.5

1 θθ

–0.5 –1

1

0

σ

σ

rr

0

1.5 r

–1 1

2

θz

0 –1 1

2

1.5 r

2

1.5 r

2

1 0

σ

σzz

1

1.5 r

1.5 r

–1 1

2

1

1

0

σ

σ

rz

rθ

0.5 0 –0.5 1

1.5 r

2

–1 1

Fig. 8. Verification for displacement boundary problem of a laminated, hollow, PZT4 cylinder. The results are the distributions of stress components (rrr , rhh , rzz , rhz , rrz , and rrh ) over the thickness at different cross-sections (z/ri = 20.0 ( ); 19.8 (); 19.6 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186 1

2183

1 0.5 0

v

u

0.5

0 –0.5 –0.5 1

1.5 r

–1 1

2

1

1.5 r

2

1.5 r

2

0

0.5 φ

w

0

–0.5

–0.5 –1 1

1.5 r

–1 1

2

Fig. 9. Verification for displacement boundary problem of a laminated, hollow, PZT4 cylinder. The results are the distributions of the displacement components (u, v, and w) and electric potential (/) over the thickness at different cross-sections (z/ri = 20.0 ( ); 19.8 (); 19.6 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



θθ

1

–0.5 –1 1

σ

σ

rr

0

1.5 r

–1 1

2

θz

0

1.5 r

1.5 r

2

1.5 r

2

rθ

1

0

σ

rz

σ

2

0 –1 1

2

1

–1 1

1.5 r

1

σ

σ

zz

1

–1 1

0

1.5 r

2

0

–1 1

Fig. 10. Verification for electric potential boundary problem of a laminated, hollow, PZT4 cylinder. The results are the distributions of stress components (rrr , rhh , rzz , rhz , rrz , and rrh ) over the thickness at different cross-sections (z/ri = 20.0 ( ); 19.8 (); 19.6 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



made of a soft elastic material, and PZT4, respectively. The material parameters for the inner layer are chosen to be E = 1.024 GPa, m = 0.3, and thus, it is approximately 100 times softer than the outer (PZT4) layer. A self-equilibrated traction rzz ¼ 9r  14 ð1 6 r 6 2Þ, rendered dimensionless with c0 ¼ 25:6 GPa, is applied at the tip-end of the

cylinder. The present results are compared with those obtained with ANSYS, because an analytical solution of this problem does not exist. The variation of the normalized non-zero stresses, displacements, and the electric potential through the thickness of the cylinder at crosssections with different axial coordinates are displayed in

2184

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

0.5

0.5

0

1 0.5 0

φ

u

w

1

–0.5

0

–0.5 –0.5 1

1.5 r

–1 1

2

1.5 r

–1

2

1

1.5 r

2

Fig. 11. Verification for electric potential boundary problem of a laminated, hollow, PZT4 cylinder. The results are the distributions of the displacement components (u and w) and electric potential (/) over the thickness at different cross-sections (z/ri = 20.0 ( ); 19.8 (); 19.6 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



0.5

0 θθ

0.5

σrr

1

σ

0

–0.5

–0.5 –1

–1

1

1.5 r

–1.5

2

1

1

1.5 r

2

1.5 r

2

0.5

0.5 σrz

σzz

0 0

–0.5

–0.5 –1

1

1.5 r

–1

2

1

Fig. 12. Verification for stress boundary problem of a hybrid laminated cylinder. The results are the distributions of non-zero stress components (rrr , rhh , rzz and rrz ) over the thickness at different cross-sections (z/ri = 0.1 ( ); 0.5 (); 1.0 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



1

0

0.5

0

–0.5

–1

1

φ

u

w

0.5

–0.5

0

1.5 r

2

–0.5

1

1.5 r

2

–1

1.5

2 r

Fig. 13. Verification for stress boundary problem of a hybrid laminated cylinder. The results are the distributions of non-zero displacement components (u and w) and electric potential (/) over the thickness at different cross-sections (z/ri = 0.1 ( ); 0.5 (); 1.0 (,)). The (dashed, dotted, etc.) lines and the symbols denote the present and ANSYS results, respectively.



Figs. 12 and 13, as computed via the two methods. Again a good agreement between the two sets of results is observed.

The real part of the lowest eigenvalue, denoted by k*, is considered for investigating the decay rate of the solution

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

as it defines the longest decay length. This value obtained for this eigenvalue is k* = 1.5447 for the hybrid laminated cylinder, whereas its value becomes 1.5452 for when the inner (softer) elastic layer is completely removed. Therefore, it may be concluded that the PZT material governs the decay rate of the solution for this hybrid laminated cylinder. When the Young’s modulus of the inner elastic material is increased 10 and 100 times, the k* reach the values of 1.4042 and 1.2507, respectively. As such, the inner elastic material has more influence on decay rate and the stresses tend to penetrate further into the interior of the cylinder when the inner layer becomes stiffer. 6. Characteristic decay length Following Miller and Horgan [20], we define the characteristic decay length L* as the length over which the stresses decay to 1% of their values at z = 0. To wit, L ¼

ln 100 ; k

ð42Þ

where k* is defined in Section 5.7. In order to illustrate the Saint-Venant end effects in piezoelectric laminated cylinders, we compute the characteristic decay lengths for two-layer, equal-thickness, hollow PZT4 cylinders with different total thicknesses and crystal orientation angles. Fig. 14 displays the computed normalized characteristic pffiffiffiffiffiffiffi decay lengths, which are defined as l ¼ L = rm t with rm ¼ ðri þ ro Þ=2 as the mean radius, for a variety of total thicknesses t ¼ ½0:25; 0:5; 0:75; 1:0; 1:25; 1:5, and crystal       orientation angles h ¼ ½0 ; 10 ; 20 ; 30 ; 40 ; 50 ;      60 ; 70 ; 80 ; 90 . It can be seen from this figure that the normalized characteristic decay length increases as the total thickness of the cylinder decreases. Fig. 14 also shows that cylinders with laminae that have the largest con trast in crystal orientation angle (i.e., h ¼ 90 ) yield the largest decay length, whereas homogeneous cylinders (i.e.,

2185



h ¼ 0 ) yield the smallest for any given total thickness. Such data can be used in designing durable and robust sensor/actuators, as the end effects and their decay length are direct metrics of delamination stresses and of deviation from uniform electromechanical fields. 7. Conclusions This study provides a semi-analytical method to analyze and quantify the end effects for laminated piezoelectric cylinders. The effects of any arbitrary self-equilibrated displacements and voltage as well as stresses and electric displacements applied at the cylinder’s ends are represented by the eigendata extracted from an algebraic eigensystem. The eigenproblem is constructed via exponentially varying displacement and voltages, that represent the solution to the homogeneous equations of equilibrium. The eigendata are used for obtaining solutions through either a leastsquares approach or the bi-orthogonality relations of the linearized eigenproblem. Comparison of present results with known analytical solutions and with three-dimensional finite element results affirms the validity of the present method. Unlike the three-dimensional finite element methods, the present approach yields direct and analytical measures of decay lengths for end effects. Furthermore, due to the reduced size of the numerical problem, the present method is more amenable for optimum cross-section design problems (e.g., minimization of delamination stresses) than three-dimensional finite element methods. An illustration of the method’s capabilities is provided in the form of analysis of characteristic decay lengths for two-layered cylinders with different crystal orientation angles and total thickness values. Acknowledgements The authors extend their gratitude to Professor Stanley B. Dong for many helpful discussions. The funding for this research was provided in part by the Academic Senate Council on Research/Faculty Grants Program (COR/ FGP) at UCLA.

4

References

l

*

3.5

3

2.5 90 0.25

60

±θ

0.5 0.75

30

1 0

1.25 1.5

t

Fig. 14. Normalized characteristic decay length (l*) of laminated PZT4 cylinders with different thicknesses (t) and crystal orientation angles (h).

[1] ANSYS, Coupled-Field Analysis Guide, third ed., ANSYS Release 5.5, September 1998. [2] H. Bai, E. Taciroglu, S.B. Dong, A.H. Shah, Elastodynamic Green’s functions for a laminated piezoelectric cylinder, Int. J. Solids Struct. 41 (2004) 6335–6350. [3] R.C. Batra, J.S. Yang, Saint-Venant’s principle in linear piezoelectricity, J. Elasticity 38 (1995) 209–218. [4] D.A. Berlincourt, D.R. Curran, H. Jaffe, Piezoelectric and piezomagnetic materials and their function in transducers, in: W.P. Mason (Ed.), Physical Acoustics: Principles and Methods, Academic Press, New York, 1964, pp. 169–270. [5] A. Borrelli, C.O. Horgan, M.C. Patria, Saint-Venant end effects in anti-plane shear for classes of linear piezoelectric materials, J. Elasticity 64 (2001) 217–236.

2186

C.W. Liu, E. Taciroglu / Comput. Methods Appl. Mech. Engrg. 196 (2007) 2173–2186

[6] A. Borrelli, C.O. Horgan, M.C. Patria, Saint-Venant’s principle for antiplane shear deformations of linear piezoelectric materials, SIAM J. Appl. Math. 62 (2002) 2027–2044. [7] A. Borrelli, C.O. Horgan, M.C. Patria, End effects for pre-stressed and pre-polarized piezoelectric solids in anti-plane shear, J. Appl. Math. Phys. (ZAMP) 54 (2003) 797–806. [8] A. Borrelli, C.O. Horgan, M.C. Patria, Exponential decay of end effects in anti-plane shear for functionally graded piezoelectric materials, Proc. Roy. Soc. Lond., Ser. A 460 (2004) 1193–1212. [9] A. Borrelli, C.O. Horgan, M.C. Patria, Saint-Venant end effects for plane deformations of linear piezoelectric solids, Int. J. Solids Struct. 43 (2005) 943–956. [10] S.B. Dong, D.B. Goetschel, Finite element analysis of edge effects in laminated composite plates, ASME J. Appl. Mech. 49 (1980) 129–135. [11] H. Fan, Decay rates in a piezoelectric strip, Int. J. Engrg. Sci. 33 (1995) 1095–1103. [12] D.B. Goetschel, T.H. Hu, Quantification of Saint-Venant’s principle for a general prismatic member, Comput. Struct. 21 (1985) 869–874. [13] C.O. Horgan, Recent developments concerning Saint-Venant’s principle: an update, Appl. Mech. Rev. 42 (1989) 295–303. [14] C.O. Horgan, Recent developments concerning Saint-Venant’s principle: A second update, Appl. Mech. Rev. 49 (1996) S101–S111. [15] C.O. Horgan, J.K. Knowles, Recent developments concerning Saint-Venant’s principle, Adv. Appl. Mech. 23 (1983) 179–269. [16] M.W. Johnson, R.W. Little, The semi-infinite strip, Q. Appl. Math. 22 (1965) 335–344. [17] M. Kazic, S.B. Dong, Analysis of restrained torsion, J. Engrg. Mech. Div. 116 (1990) 870–891.

[18] J.K. Knowles, On Saint-Venant’s principle in the two-dimensional linear theory of elasticity, Arch. Ration. Mech. Anal. 21 (1966) 1– 22. [19] R.W. Little, S.B. Childs, Elastostatic boundary region problem in solid cylinders, Q. Appl. Math. 25 (1967) 261–274. [20] K.L. Miller, C.O. Horgan, Saint-Venant end effects for plane deformations of elastic composites, Mech. Compos. Mater. Struct. 2 (1995) 203–214. [21] X. Ruan, S.C. Danforth, A. Safari, T.W. Chou, Saint-Venant end effects in piezoceramic materials, Int. J. Solids Struct. 37 (2000) 2625– 2637. [22] E. Taciroglu, C.W. Liu, S.B. Dong, C.K. Chun, Analysis of laminated piezoelectric circular cylinders under axisymmetric mechanical and electrical loads with a semi-analytic finite element method, Int. J. Solids Struct. 41 (2004) 5185–5208. [23] E. Taciroglu, C.W. Liu, Analysis and design of multi-modal piezoelectric layered tubular sensors and actuators, Smart Mater. Struct. 14 (2005) 605–614. [24] J.Q. Tarn, L.J. Huang, Saint-Venant end effects in multilayered piezoelectric laminates, Int. J. Solids Struct. 39 (2002) 4979–4998. [25] R.A. Toupin, Saint-Venant’s principle, Arch. Ration. Mech. Anal. 18 (1965) 83–96. [26] S. Vel, R.C. Batra, Analytical solutions of rectangular thick laminated plates subjected to arbitrary boundary conditions, AIAA J. 37 (1999) 1464–1473. [27] S. Vel, R.C. Batra, The generalized plane strain deformations of thick anisotropic composite laminated plates, Int. J. Solids Struct. 37 (2000) 715–733.