Boundary Element Modeling of Dynamic Bending of a Circular Piezoelectric Plate

Accurate modelling of a coupled dynamic electro-mechanical response of circular piezoelectric plates under various loading conditions is of particular importance. Piezoelectric plates are not only basic structural elements, but with certain considerations can be conveniently fit for numerical simulation of piezoelectric sensors and transducers. In this work, a Laplace domain direct boundary element formulation is applied for dynamic analysis of three-dimensional linear piezoelectric moderately thick circular plates. Zero initial conditions, vanishing body forces and the absence of the free electrical charges are assumed. Weakly singular expressions of Laplace domain boundary integral equations for the generalized displacements are employed. Spatial discretization is based on the nodal collocation method. Mixed boundary elements are implemented. The geometry of the elements, generalized displacement and generalized tractions are represented with different shape functions: quadratic, linear and constant, accordingly. Integral expressions of the three-dimensional Laplace domain piezoelectric displacement fundamental solutions are used. After solving the problem on a set of Laplace transform parameter values, time-domain solutions are retrieved from the corresponding Laplace domain solutions by employing a numerical inversion routine. Numerical example is provided to show reliability and accuracy of the proposed boundary element formulation.


Introduction
Due to the direct piezoelectric and converse piezoelectric effects, piezoelectric materials have a very practically significant ability to convert the mechanical energy into the electric energy and vice versa. Circular piezoelectric plates are extensively used in microelectromechanical devices, transducers, sensors, actuators, transformers, headphones, loudspeakers, etc. [1][2][3]. Analytical approaches for the time-harmonic and especially transient analysis of piezoelectric structures are limited to simple boundary conditions and specific geometries [4][5][6].
Therefore, suitable numerical techniques are required to predict the dynamic behaviour of the general piezoelectric structures, to improve their efficiency, reliability and service life.
The numerical solution of the three-dimensional piezoelectric static, time-harmonic and transient problems is usually carried out using the Finite Element Method (FEM). The Boundary Element Method (BEM) is considered to be a good alternative to FEM for problems of small deformations and with linear constitutive material behaviour. Typically, BEM requires to discretize only the boundary of the considered region which is advantageous compared to the FEM. Though various BEM formulations has been successfully applied for linear isotropic elastodynamic problems [7][8][9][10][11][12], number of works dedicated to application of BEM to the three-dimensional anisotropic linear elastic and piezoelectric dynamic problems is still low. This is due to the fact that the usability of the conventional boundary element formulation depends on the availability of the fundamental solutions. In the case of the linear piezoelectricity, as well as the linear anisotropic elasticity, the static and dynamic fundamental solutions are not available in the explicit closed-form expressions. To avoid the complications, associated with a numerical implementation of the dynamic fundamental solutions, a Dual Reciprocity Boundary Element Method (DRBEM) which employs only the static fundamental solutions has been developed [13]. However, to correctly take into account the inertial effects in DRBEM solution procedure, additional collocation nodes are required to be placed inside the considered domain.
In this paper, a transient dynamic behaviour of three-dimensional linear piezoelectric moderately thick circular plates is modelled by the Laplace domain direct boundary element formulation. Zero body forces, zero initial conditions and the absence of the free electrical charges are assumed. Weakly singular representation of Laplace domain boundary integral equations for the generalized displacements is used. Nodal collocation method is used for the spatial discretization. To approximate the geometry of the elements, generalized displacement and generalized tractions the idea of the mixed boundary elements is employed: geometry is interpolated with the quadratic shape functions, displacements and tractions are interpolated with the linear and constant shape functions, respectively. Three-dimensional Laplace domain piezoelectric

Basic equations
We consider a three-dimensional, homogeneous, linear piezoelectric solid with volume 3 R Ω ⊂ and boundary . Γ = ∂Ω Under quasi-electrostatic assumption, without free electrical charges, with zero body forces and zero initial conditions, the dynamic problem of a linear piezoelectric solid in the Laplace domain is described by following coupled differential equations written in terms of generalized quantities: where , 1, 3, i l = , 1, 4; j k = ρ and s represent the mass density and the complex-valued Laplace transform parameter; k U and ijkl C the generalized displacements and generalized elasticity tensor, which are defined by using contracted notation [14] as follows , , , , 1, 3, , , , 1, 3, 4, , , , 1, 3, 4, , , 1, 3, 4.
In the Eqs. Taking into account the symmetries of these tensors the following symmetry relation for the generalized elasticity tensor can be established: .
ijkl lkji The boundary conditions are expressed as with the generalized tractions being defined as where j t and q denote the elastic tractions and normal charge flux density; jk σ and k D stand for the Cauchy stress tensor and electric displacements; k n is the components of the outward unit normal vector to the boundary. U Γ is the part of the external boundary Γ where the generalized displacements k U * are prescribed and T Γ is the part of the boundary where the generalized tractions k T * are given.

Boundary element formulation
In order to apply the boundary element method, boundary integral equations (BIEs) for the boundary value problem in the Laplace domain defined in Eqs.
(1)-(2) and Eqs. (6)- (7) is formulated. In the present work, the weakly-singular representation of the BIEs for the generalized displacements k U is employed: where y is the field point and x is the source point. The boundary Γ is then discretized into the finite number of quadrangular boundary elements. The mixed boundary elements are used for the spatial discretization of the boundary integral Eq. (9) which means that the geometry of the elements and the variation of the generalized displacements and generalized tractions over each element are interpolated by different shape functions: quadratic, linear and constant, respectively. By carrying out the standard nodal collocation procedure, necessary integration and applying the known boundary conditions, a complex-valued system of linear algebraic equations is obtained for a selected value of the Laplace transform parameter s:  (12) , = − r y x , r = r (13) where the superscripts "S" and "R" respectively indicate the singular (static) and regular (dynamic) parts of the corresponding generalized piezoelectric fundamental solutions.
The generalized traction fundamental solution is related to the generalized displacement fundamental solution as follows where ( ) i n y is the outward unit normal vector to the boundary at the field point y.
Dynamic part of the Laplace domain piezoelectric three-dimensional generalized displacement fundamental solution can be expressed by a two-dimensional integral over the half of a unit sphere as [15]  n (15) while the static part is expressed as the integral over a unit circle as

Laplace transform inversion
In present formulation we use the convolution quadrature method (CQM) [16,17] based on the Radau IIA Runge-Kutta method (RKCQM) [18,19] reformulated as a numerical method for inverting a Laplace transform. Assuming that RKCQM is given by where N is the total number of the time steps, t ∆ is the time step size and 0 1 R < < is a CQM parameter.
A circular piezoelectric plate with the diameter 0.01 m d = and thickness 0.00125 m h = is considered ( Fig. 2). At its upper surface 3 2 x h = + the plate is subjected (Fig. 3) to a uniform normal traction load 3 ( ), t t H t * = 9 10 Pa t * = − with ( ) H t being the Heaviside function. The rest of the plate's surface is free of the generalized tractions. The plate is clamped on its circumferential edge and electrically grounded on the whole surface. We consider the PVDF piezoelectric material with mass density  h ) m. It can be clearly observed that the results on the most coarse mesh and most refined mesh differ only slightly. The comparison of the boundary element and finite element solutions (Fig. 5)

Concluding remarks
The Laplace domain direct boundary element formulation for dynamic analysis of three-dimensional linear piezoelectric moderately thick circular plates is presented in this paper. Dynamic problems of linear piezoelectric solids are considered under quasi-electrostatic assumption, without body forces, in the absence of free electrical charges and with vanishing initial conditions. BEM formulation is based on the weakly-singular representation of the BIEs for the generalized displacements. The spatial discretization is done with the mixed interpolation of the geometry of the elements and the boundary fields. Dynamic part of the Laplace domain piezoelectric three-dimensional generalized displacement fundamental solution is expressed by a two-dimensional integral over the half of a unit sphere and the static part is expressed as the integral over a unit circle. Convolution quadrature method based on the Radau IIA Runge-Kutta method is employed to recover time-domain solutions from the Laplace domain solutions. Numerical example is provided to verify the validity of the proposed boundary element approach. Comparison of the present results with those from the finite element analysis confirmed accuracy of the formulation.