Abstract
In this study, a piezoelectric laminate is analyzed by applying the Laplace transform and its numerical inversion, Fourier transform, differential quadrature method (DQM), and state space method. Based on the modified variation principle for the piezoelectric laminate, the fundamental equations for dynamic problems are derived. The solutions for the displacement, stress, electric potential, and dielectric displacement are obtained using the proposed method. Durbin’s inversion method for the Laplace transform is adopted. Four boundary conditions are discussed through the DQM. The proposed method is validated by comparing the results with those of the finite element method (FEM). Moreover, this semianalytical method is further extended to describe the dynamic response of piezoelectric laminated plates subjected to underwater shocks by introducing Taylor’s fluid-structure interaction algorithm. Both air-backed and water-backed laminated plates are investigated, and the effect of the fluid is examined. In the time domain, the electric potential and displacements of sample points are calculated under four boundary conditions. The present method is shown to be accurate and can be a useful method to calculate the dynamic response of underwater laminated sensors.
1. Introduction
Piezoelectric materials are extensively used as sensing elements in various engineering fields, such as aerospace engineering, ocean engineering, civil engineering, and mechanical engineering. In addition, piezoelectric materials are one of the most important components in smart structures which can sense and drive changes in external circumstances and then be modified in response by computers and processors [1]. Since the Curie brothers first discovered the direct piezoelectric effect in single-crystal quartz in 1880 and Gabriel Lippmann discovered the converse piezoelectric effect in 1881, various types of piezoelectric materials have been found, including both natural materials and composite materials [2, 3].
Composite laminates are a common smart structure. As shown in Figure 1, the structure is composed of stacked layers of different materials; the individual layers are generally orthotropic or transversely isotropic. Owing to the various characteristics of the individual layers, a hybrid laminate may exhibit multiple functions. Piezoelectric laminates are usually piled as sensors and actuators that will react to or induce a displacement or electric potential [4, 5].
Pioneering research on piezoelectric laminates has been conducted primarily using finite element models (FEMs) and the space state method. A FEM has stable convergence and its solution is effective; however, it requires complex grids and long computation times to produce results with high accuracy [6, 7]. Mukherjee and Saha [8] focused on large deformations of piezoelectric beams using a FEM based on the first-order shear deformation theory and the Newton–Raphson method. Garcia et al. [9] developed a layer-wise FEM for a piezoelectric plate using the Reissner mixed variational principle. Dash and Singh [10] employed a C0 isoparametric FEM to solve the nonlinear free vibration problem. Wankhade and Bajoria [11] considered the higher-order shear deformation theory to carry out a FEM of a piezoelectric beam under static and dynamic excitations.
On the other hand, the state space method is also an efficient and accurate method. In view of the relationships among transverse vectors, it is suitable for laminated structures [12]. Different approaches have been investigated in the state space method for piezoelectric laminates. D’Ottavio and Kröplin [13] combined the Reissner variational statement and constitutive equations to simulate piezoelectric laminates. Tiersten [14] demonstrated Hamilton’s principle for linear piezoelectric media under different boundary conditions and discussed both homogeneous and inhomogeneous problems. In addition, Qing et al. [15] proposed a modified mixed Hellinger–Reissner (H–R) variational principle based on the work of Steele and Kim [16].
The state space method is the primary focus of this study. The variational principle for the state space method is a complex differential equation including several variables, and thus some transformations should be employed. The Laplace transform and differential quadrature method (DQM) are the two main transformations used in this study.
The Laplace transform is usually used to change a function in the time domain to one in a complex frequency field for dynamic problems [17]. For this reason, in recent studies, the Laplace transform has generally been adopted to analyze both statics and dynamics problems, and the numerical inversion of the Laplace transform has been studied. Zhao [18] proposed two algorithms that divide the integration interval into small subspaces with different integration steps. Durbin [19] combined a trigonometric series with an inversion series and presented rather smaller results. Qing et al. [15] subsequently made comparisons to the inverse algorithms proposed by both Zhao and Durbin and improved Zhao’s algorithms by using Subbotin splines. Wang et al. [6] improved the calculation efficiency with two novel algorithms and additionally compensated for the deviation node in Qing’s algorithms. Durbin’s inversion method is used in this study owing to its high efficiency, and the results are sufficiently precise.
First proposed in 1971 by Bellman et al. [20] based on integral quadrature, DQM is another method for dealing with partial differential factors. Moreover, with DQM, complex boundary conditions can be analyzed. Extensive research has focused on the choice of grid points, boundary conditions, weight coefficients, and the convergence behavior. Several researchers [21–23] have eliminated the unreliable results from uniform grid points and offered various distribution forms that could satisfy higher orders and different simulations. Chen et al. [24] employed the Chebyshev polynomial to obtain the weighting coefficients of the matrix form. To date, DQM has been widely applied for composite structures. Du and Shu [23] chose DQM to deal with the clamped and simply supported boundary conditions of beams and plates. Pradhan and Murmu [25] used DQM to analyze the vibration of functionally graded material beams.
Piezoelectric laminates are employed for numerous applications in the sea. It is important for underwater structures to consider the interaction between the water and the structure. There are few studies regarding underwater piezoelectric laminates, and thus the fluid-structure interaction (FSI) is also considered in this study. Taylor [26] proposed a one-dimensional model to analyze the FSI. It has been further verified under various conditions with different materials, fluid media, and structures. Schiffer and Tagarielli [27] examined the FSI response with various materials and structural combinations. Blom [28] investigated a piston problem and combined the one-dimensional equation and the Euler equation to propose an algorithm. Deshpande and Fleck [29, 30] studied sandwich beams subjected to underwater waves using a lumped parameter model and FEM.
In this study, we discuss three-dimensional semianalytical solutions for piezoelectric laminates subjected to an underwater shock. The state space method is established based on the constitutive equations for piezoelectric materials. Both air-backed and water-backed laminated plates are investigated based on the effect of the FSI. Combined with the Fourier transform, DQM, and Laplace transform, the present model is analyzed under four different boundary conditions. The results are verified using finite element analysis software.
2. Mathematical Description
2.1. Fundamental Equations
A piezoelectric laminate consists of a series of layers of different materials, including piezoelectric materials. Consider a laminated plate comprising three layers of PZT/Gr-Epoxy/PZT materials, where PZT means lead zirconate titanate piezoelectric ceramics. As shown in Figure 1, the laminate is subject to an underwater shock. Thus, the upper layer is always exposed to the water, while the lower layer may be exposed to either water or air.
Briefly, the piezoelectric effect is demonstrated to change polarization when a mechanical stress is applied or to create a mechanical deformation with the introduction of additional voltage. In view of basic theory, the constitutive equation for piezoelectric materials can be written simply as the following form:where , , , and are the stress, strain, electric field, and electric displacement, respectively, and , , and are the elastic, piezoelectric, and dielectric coefficients, respectively. The subscripts p, q, i, and j denote different directions in the material coordinates (i, j = 1, 2, 3; p, q = 1, 2, 3, 4, 5, 6). The superscript E indicates parameters measured at a constant electric field, and the superscript S denotes those measured at a constant or zero strain field. The superscripts will be omitted in the subsequent derivation process.
Assuming that the piezoelectric materials are orthotropic or have orthotropic symmetry relative to the x-y coordinate plane, constitutive equation (1) can be extended to the following matrix form:in the Cartesian coordinate system; p, q = x, y, z, yz, xz, xy, and i, j = x, y, z. For the model of a piezoelectric laminate shown in Figure 1, the subsequent equations and derivations will be described in the Cartesian coordinate system.
In terms of geometric relationships, the relationships between the strains and displacements can be described in the form of gradient equations:where u, , and denote the displacements in the three directions of the Cartesian axes; and α, β, and γ are partial differential operators of the three Cartesian axes with respect to the following variable: .
In light of the quasi-static Maxwell equations, the electric fields, , and electric potential, , have similar relationships as the strains and displacements:
Two general designations are defined as and for the two groups including four dual vectors, while the remaining vectors of the variables are given as , , and . Thus, the subsequent derivations can be simplified with these five designations. Constitutive equation (2) can then be rewritten as follows:where , , , and represent the transformations of the 9 × 9 matrix in equation (2) and .
Hence, equations (3) and (4) become the two following matrix forms:where .
According to the Reissner mixed variational theorem, Reissner’s energy density function is given by the above matrices:
Owing to the special form of orthotropic constitutive equations, the secondary vectors in equation (7) are replaced, and the integral parts are derived to dot-product forms:where for dynamic problems in the time domain.
In terms of the Reissner variational principle, the energy equation with an internal force is given by the following equation:where h indicates the thickness of the layer and is the area of the layer.
Therefore, based on equations (8) and (9), the state space method form can be deduced as follows:where , , , and .
In addition, the in-plane vectors and rest vectors can be evaluated by the expression in the following equation:
The partial differential operator of z in equation (10) can be eliminated by transforming the coefficient matrix into an exponential matrix as follows:where k indicates the kth layer and represents the value of the thickness at the kth layer.
Therefore, the state vectors of the piezoelectric laminate at z can be written as follows:where
2.2. Fluid-Structure Interaction (FSI)
As shown in Figure 1, the piezoelectric laminate is subject to an underwater explosion, resulting in shock waves on both the upper and lower surfaces of the laminate. In accordance with the one-dimensional FSI model proposed by Taylor [26] in 1963, the effect of the shock waves on the laminate can be calculated.
The shock wave is simply described as a function of the time, t, and the vertical distance to the plane, z, which in this paper is equal to zero:where denotes the wave velocity in water as a constant and and are empirical formulas for the explosive charge weight, (kg), and perpendicular distance from the upper surface of the laminate, S (m), respectively. TNT is used for the explosion in this study:
Under these conditions, the planes influenced by the shock wave include both the upper and lower surfaces of the laminate [29, 30]. The wave pressure on one plane consists mainly of the primary shock wave , reflected wave, , and rarefaction wave, :where is the density of water and .
Consequently, when the laminate is completely submerged underwater, the wave pressures on the upper and lower surfaces are given by the following:
In this paper, the interactions between the fluid and structure are divided into two conditions. As the wave pressure on the upper surface is caused directly by the underwater shock, the lower surface is exposed to both water and air, which here are referred to as water-backed and air-backed laminated plates. Because the density of air is far smaller than the density of water, the wave pressure on the air-backed surface is assumed to be zero:
Hence, on the upper and lower surfaces of the laminate, the boundary conditions are confirmed as follows according to the FSI:where for the water-backed plate and for the air-backed plate.
3. Boundary Conditions and Transform Methods
3.1. Laplace Transform Methods
The Laplace transform method is adopted to obtain the response in the time domain, which is described in the two following equations:
In this paper, Durbin’s inversion method given by equation (23) is adopted. The trigonometric series is obtained, and the error of the results is independent of t, and thus the analytical model can yield good results with high efficiency:where , , and is the real part of imaginary number , and , where denotes the duration of the model in the time domain. Here, the value of K is 100.
3.2. Fourier Transform Methods
To represent different means of support, the boundary conditions at x = 0, a are expressed in four ways, as given by equations (24)–(32).
Simply supported (x = 0)-simply supported (x = a) (S-S):
Clamped (x = 0)-clamped (x = a) (C-C):
Clamped (x = 0)-free (x = a) (C-F):
Clamped (x = 0)-simply supported (x = a) (C-S):
The Fourier transform method is employed here to change the form of . According to the coefficient matrix in equation (10), the Fourier transform indicates the following:
In addition, the Fourier transform of the pressure on the surfaces yieldswhere
3.3. Differential Quadrature Method (DQM)
While the Fourier transform method is used to change the form of , DQM is used to parse the form of .
By selecting some grid points in the domain of x, the differential form of a continuous function, , can be described by an approximate function that is calculated with the values of all grid points and a weighting coefficient matrix [7, 31, 32]:where the superscript i indicates the ith-order derivative, M is the total number of sampling points, and denotes the weighting coefficients.
The Chebyshev–Gauss–Lobatto points are chosen for their efficient speed of convergence and high computational accuracy [7]. Therefore, the sample points of the x-coordinate can be chosen as follows:
3.4. Normalization of the Variables
The variables need to be normalized to prevent nonconvergence of the computation caused by large gaps in magnitude between different parameters.
According to the constitutive equation in equation (1), the variables can be normalized in the following ways:where .
For FSI, the density of water, , and wave velocity underwater, , are normalized as follows:
The of the normalized variables will be omitted in the following equations.
In view of the semianalytical transformations mentioned above, for the boundary conditions of C-C, S-S, C-F, and C-S, the state space method form in equation (10) is further deduced as follows:where the subscript i indicates the ith grid point in the DQM (i = 1, 2, …, M), , , and .
Due to DQM, the vectors of the four support means as equations (26)–(33) can be rewritten as
Based on equation (11) and the boundary conditions as stated above, equation (39) should be deduced to avoid singular matrices in it, which means the coefficient matrix in equation (10) should be deduced into different forms for the four boundary conditions:where m, n = 2∼M − 1 in the matrix andwhere m = 2∼M − 1, n = 1∼M − 1 in the matrix, andwhere m, n = 1∼M in the matrix andwherewhere
4. Numerical Results
The model in Figure 1 has two components: the first and third layers are composed of PZT, and the second layer is Gr/Epoxy. All the material parameters are listed in Table 1, including the elastic modulus, density, piezoelectric constant, and dielectric constant. The length, width, and thickness of the model are represented as a, b, and h, respectively, and each layer has a thickness of h/3. The values of a, b, and h are provided in Table 2 (see Tables 3 and 4).
The semianalytical underwater laminated model based on the state space method incorporating the Laplace transform, Fourier transform, and DQM is programmed using the Mathematica platform.
4.1. Comparison of the Present Method and FEM without FSI
First, the piezoelectric laminate under mechanical pressure is analyzed without considering the FSI. A sinusoidal transverse load, q0, is applied to the upper surface of the model, where . For the electrical boundary, the electric potential, , is set to zero at both the upper and lower surfaces. Meanwhile, the numerical results are compared with the simulation with the FEM.
Taking the C-C boundary condition as an example, the results of both methods are shown in Figures 2(a) and 2(b); the chosen sample points are x = a/2, y = b/2 in the x-y plane along the thickness, z. Figure 3 shows the distribution of the electric potential in the x-z plane at y = b/2. Figures 2 and 3 represent the results at t = 0.5 s, while Figure 4 displays the time-history curves of the vectors of the sample point at x = a/2, y = b/2, and z = h/2.
(a)
(b)
(a)
(b)
(a)
(b)
As shown in Figure 2, compared with the FEM results, the values of (m) and ϕ (V) of the sample points obtained with the present method have similar trends along the thickness; the correlation coefficient in Figure 2(a) is 99.997%, while that in Figure 2(b) is 99.999%. The correlation coefficient can be calculated as in the following equation:
Figure 3 shows that both the FEM and present method produce similar distributions, and Figure 4 shows almost the same trends over time obtained with the two methods. These similar trends, distributions, and correlation coefficients confirm the reliability of the present method.
4.2. Comparison of the Present Method and FEM with FSI
The dynamic responses of underwater piezoelectric laminated plates are verified by FEM. Taking the air-backed C-C boundary condition into consideration, the upper layer of the plate is exposed to the water. Therefore, a spring foundation is set on the upper layer in the FEM model, where the spring constant is zero, and the viscous damping here is , which shows the underwater properties. Equation (48) denotes the relationship of the spring constant and viscous damping. Besides, the time step is set as 0.0001 s due to equation (16) of shock wave:where is spring constant and is viscous damping.
The same sample point at x = a/2, y = b/2, and z = h/2 in the time domain is chosen. Figures 5 and 6 show the comparisons between two curves of both present method and FEM. The chosen period of FEM curves is from 0 s to 0.3 s.
(a)
(b)
(a)
(b)
It is clear that the time history of and the peak values of ϕ are close, but the curves of ϕ show a sharp fluctuation during the first 0.03 s in both the present method and FEM. The FEM takes more steps at first 0.03 seconds to analyze the shock wave; besides the values are more unstable. While the present method owns a stable change during the first 0.003 s, the values show a slight fluctuation in the time domain due to the limitation of computational accuracy.
4.3. Influence of the FSI Based on the Present Method
The FSI is coupled to the piezoelectric laminate in this section. The vectors on the upper and lower surfaces are applied as equations (19) and (20). In addition, the plane stress on the surfaces is defined as zero to satisfy the one-dimensional FSI theory. To investigate the influence of the FSI, air-backed and water-backed laminated plates are simulated for each boundary condition. The electric vectors are defined as ϕ = 0 on the upper surface and the lower surface.
On account of the equations for the shock wave given in equations (15) and (16), for and , = 1 kg and S = 500 m to establish the boundary conditions on the upper and lower surfaces. The results for the four boundary conditions of S-S, C-C, C-F, and C-S are shown in Figures 7–14. In Figures 7–10, t = 0.05 s, and sample points of x = a/2 and y = b/2 along the thickness, z, are chosen. The values on the boundaries of some figures are shown in white, which indicates that they are much larger or much smaller than the rest of the figure.
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
For various boundary conditions, both the distributions and values of the electric potential at the x-z plane are different. For the C-C, C-F, and C-S boundary conditions, the distributions of the values are almost the same under both air-backed and water-backed conditions, and the FSI on the surfaces lessens the external force, which leads to a reduction in ϕ. For the S-S condition, the distribution of value shows an obvious difference; the maximum values shown in the figures are at different positions, which is because there are rotational degrees of freedom and the stress tends to disperse at 1/4 and 3/4 of the x-axis under the smaller external force. The same phenomenon can also be observed in the C-S boundary condition, where the large values disperse near the simply supported boundary. For the C-F condition, the stress concentrates on the clamped boundary and the values are obviously large, while the values close to the free boundary are almost 0.
The response time history of each vector can also be calculated. Here, the sample point at x = a/2, y = b/2, and z = h/2 in the time domain is chosen for analysis.
In Figures 11–14, it is clear that the values of each vector at the chosen point are lessened by the rarefaction wave on the lower surface, and the values under water-backed condition are about half of the values under air-backed condition. The decreasing of water-backed value seems to be slower than that of air-backed value. The water-backed condition shows more buffering against the shock wave. The four kinds of boundary conditions indicate the close value of , and ϕ under S-S condition seems to be the smallest due to the more degrees of freedom on the boundary.
5. Discussion and Conclusions
A novel and efficient method is proposed in this paper to solve the problem of the dynamic responses of underwater piezoelectric laminated plates. This method first combines the Laplace transform, Fourier transform, and DQM to change the differential forms of x, y, and t. The state space method is then introduced to deduce the fundamental equations of the piezoelectric material. Finally, the FSI is employed to couple the piezoelectric laminate and influence of the water.
In this study, the finite element method is chosen to verify the results of the present method for the piezoelectric problem. Although there is some deviation in the values of some vectors, the correlation coefficients of the results and overall trends indicate that the present method is accurate and reliable. For the comparison of FSI on the C-C condition, the FEM dealing with shock load reveals instability, and the present method has its advantage in the dynamic response.
As the FEM has limits in coupling the piezoelectric-mechanical-fluid problem, the one-dimensional FSI is used to analyze this condition. The results show that the present method has good convergence in conveying external forces such as a shock wave in the form of the FSI. Additionally, the FSI has an obvious buffer effect on the external force.
When the underwater shock was applied on the model under different boundary conditions, it can be easily found that the water on the surfaces has buffering to the flow. Different boundary conditions show their specific distributions of the vectors, which can be the basis of sensors under different supporting boundaries.
The method proposed in this paper aims to analyze underwater piezoelectric laminates, which can be applied as sensors and actuators in marine engineering. Furthermore, it also provides a new semianalytical solution to multiphysics problems.
Data Availability
The material and model data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The authors gratefully acknowledge financial support from the National Natural Science Foundation of China (Grant nos. 51879231, 51679214, and 51979244).