Analysis and active control of geometrically nonlinear responses of smart FG porous plates with graphene nanoplatelets reinforcement based on B\'ezier extraction of NURBS

TL;DR

This paper introduces a computational method combining isogeometric analysis and higher-order shear deformation theory to analyze and actively control the nonlinear responses of functionally graded porous plates reinforced with graphene nanoplatelets and integrated with piezoelectric layers.

Contribution

It develops a novel Bézier extraction-based isogeometric analysis approach for nonlinear analysis and active control of FG porous plates with GPLs and piezoelectric layers.

Findings

The method accurately predicts nonlinear behaviors under various conditions.

Active control effectively suppresses nonlinear responses.

Porosity and GPL distribution significantly influence structural responses.

Abstract

In this paper, we propose an effective computational approach to analyze and active control of geometrically nonlinear responses of functionally graded (FG) porous plates with graphene nanoplatelets (GPLs) reinforcement integrated with piezoelectric layers. The key concept behind this work is to utilize isogeometric analysis (IGA) based on B\'ezier extraction technique and $C^0$-type higher-order shear deformation theory ($C^0$-HSDT). By applying B\'ezier extraction, the original Non-Uniform Rational B-Spline (NURBS) control meshes can be transformed into B\'ezier elements which allow us to inherit the standard numerical procedure like the standard finite element method (FEM). In this scenario, the approximation of mechanical displacement field is calculated via $C^0$-HSDT whilst the electric potential field is considered as a linear function across the thickness of each piezoelectric sublayer. The FG plate includes internal pores and GPLs dispersed into metal matrix either uniformly or non-uniformly along plate's thickness. To control responses of structures, the top and bottom surfaces of FG plate are firmly bonded with piezoelectric layers which are considered as sensor and actuator layers. The geometrically nonlinear equations are solved by Newton-Raphson iterative procedure and Newmark's integration. The influence of porosity coefficient, weight fraction of GPLs as well as external electrical voltage on geometrically nonlinear behaviors of plate structures with various distributions of porosity and GPLs are thoroughly investigated. A constant displacement and velocity feedback control approaches are then adopted to actively control geometrically nonlinear static and dynamic responses, where structural damping effect is taken into account, based on a closed-loop control with sensor and actuator layers.

Tables

Table 1. Table 1 : Normalized central deflection w ¯ ¯ 𝑤 \overline{w} of a CCCC isotropic square plate ( a / h = 100 𝑎 ℎ 100 a/h=100 ) under a uniform load ( P = q 0 ​ a 4 / ( E ​ h 2 ) 𝑃 subscript 𝑞 0 superscript 𝑎 4 𝐸 superscript ℎ 2 P=q_{0}a^{4}/(Eh^{2}) ).

$P$	Present	Analytical [56]	MXFEM [57]	IGA-RPT [58]
17.79	0.2348	0.237	0.2328	0.2365
38.3	0.4663	0.471	0.4738	0.4692
63.4	0.6873	0.695	0.6965	0.6908
95.0	0.8983	0.912	0.9087	0.9024
134.9	1.1016	1.121	1.1130	1.1060
184.0	1.2960	1.323	1.3080	1.3008
245.0	1.4875	1.521	1.5010	1.4926
318.0	1.6728	1.714	1.6880	1.6784
402.0	1.8492	1.902	1.8660	1.8552

Table 2. Table 2 : Material properties of the metal and piezoelectric materials.

Properties	Core layer					Piezoelectric layer
	Ti-6Al-4V	Aluminum oxide	Copper	GPLs		PZT-G1195N
Elastic properties
$E_{11}$ (GPa)	105.70	320.24	130	1010		63.0
$E_{22}$ (GPa)	105.70	320.24	130	1010		63.0
$E_{33}$ (GPa)	105.70	320.24	130	1010		63.0
$G_{12}$ (GPa)	$-$	$-$	$-$	$-$		24.2
$G_{13}$ (GPa)	$-$	$-$	$-$	$-$		24.2
$G_{23}$ (GPa)	$-$	$-$	$-$	$-$		24.2
${\nu }_{12}$	0.2981	0.26	0.34	0.186		0.30
${\nu }_{13}$	0.2981	0.26	0.34	0.186		0.30
${\nu }_{23}$	0.2981	0.26	0.34	0.186		0.30
Mass density
$\rho$ (kg/m3)	4429	3750	8960	1062.5		7600
Piezoelectric coefficients
$k_{31}$ (m/V)	$-$	$-$	$-$	$-$		$254\times {10}^{-12}$
$k_{32}$ (m/V)	$-$	$-$	$-$	$-$		$254\times {10}^{-12}$
Electric permittivity
$p_{11}$ (F/m)	$-$	$-$	$-$	$-$		$15.3\times {10}^{-9}$
$p_{22}$ (F/m)	$-$	$-$	$-$	$-$		$15.3\times {10}^{-9}$
$p_{33}$ (F/m)	$-$	$-$	$-$	$-$		$15.3\times {10}^{-9}$

Table 3. Table 3 : Tip node deflection of cantilever piezoelectric FGM (Aluminum oxide/Ti-6Al-4V) plate ( a = b = 0.4 𝑎 𝑏 0.4 a=b=0.4 m, h c = 5 subscript ℎ 𝑐 5 h_{c}=5 mm and h p = 0.1 subscript ℎ 𝑝 0.1 h_{p}=0.1 mm) under a uniform load ( q 0 = 100 subscript 𝑞 0 100 q_{0}=100 N/m 2 ) and various input voltages ( × 10 − 4 absent superscript 10 4 \times 10^{-4} m).

$n$	Method	Input voltages (V)
		0	20	40
$n=0$	Present	-2.5437	-1.3328	-0.1229
	CS-DSG3 [60]	-2.5460	-1.3346	-0.1232
$n=0.5$	Present	-1.6169	-0.8418	-0.0667
	CS-DSG3 [60]	-1.6199	-0.8440	-0.0681
$n=5$	Present	-1.1233	-0.5808	-0.0382
	CS-DSG3 [60]	-1.1266	-0.5820	-0.0375
$n=\infty$	Present	-0.8946	-0.4608	-0.0271
	CS-DSG3 [60]	-0.8947	-0.4609	-0.0271

Equations

$$\left\{\begin{array}[]{l}E\left(z\right)={E_{1}}\left[{1-{e_{0}}\lambda\left(z\right)}\right],\\ G\left(z\right)=E\left(z\right)/\left[{2\left({1+\nu\left(z\right)}\right)}\right],\\ \rho\left(z\right)={\rho_{1}}\left[{1-{e_{m}}\lambda\left(z\right)}\right],\end{array}\right.$$

$$\lambda\left(z\right)=\left\{{\begin{array}[]{*{20}{l}}{\textrm{cos}\left({\pi z/h_{c}}\right)},&\textrm{Porosity distribution 1}\\ {\textrm{cos}\left({\pi z/2h_{c}+\pi/4}\right)},&\textrm{Porosity distribution 2}\end{array}}\right.$$

$$e_0=1-\frac{E_2^{\prime }}{E_1^{\prime }},(0\le e_0<1),$$

$$\frac{E\left(z\right)}{E_1}={\left(\frac{\rho \left(z\right)/{\rho }_1+0.121}{1.121}\right)}^{2.3}\text{for}\left(0.15<\frac{\rho \left(z\right)}{{\rho }_1}<1\right).$$

$$e_m=\frac{1.121\left(1-\sqrt[2.3]{1-e_0\lambda \left(z\right)}\right)}{\lambda \left(z\right)}.$$

$$\nu \left(z\right)=0.221\varsigma +{\nu }_1\left(0.342\varsigma {}^2-1.21\varsigma +1\right),$$

$$\varsigma =1.121\left(1-\sqrt[2.3]{1-e_0\lambda \left(z\right)}\right).$$

$${V_{GPL}}=\left\{{\begin{array}[]{*{20}{l}}{{S_{i1}}\left[{1-\cos\left({\pi z/h_{c}}\right)}\right]},&\textrm{Pattern A}\\ {{S_{i2}}\left[{1-\cos\left({\pi z/2h_{c}+\pi/4}\right)}\right]},&\textrm{Pattern B}\\ {{S_{i3}}},&\textrm{Pattern C}\end{array}}\right.$$

$$\frac{{\Lambda }_{GPL}{\rho }_m}{{\Lambda }_{GPL}{\rho }_m+{\rho }_{GPL}-{\Lambda }_{GPL}{\rho }_{GPL}}\times {\int }_{-h_c/2}^{h_c/2}\left[1-e_m\lambda \left(z\right)\right]dz={\int }_{-h_c/2}^{h_c/2}{V}_{GPL}\left[1-e_m\lambda \left(z\right)\right]dz.$$

$$E_1=\frac{3}{8}\left(\frac{1+{\xi }_L{\eta }_L{V}_{GPL}}{1-{\eta }_L{V}_{GPL}}\right)E_m+\frac{5}{8}\left(\frac{1+{\xi }_W{\eta }_W{V}_{GPL}}{1-{\eta }_W{V}_{GPL}}\right)E_m,$$

$$\begin{array}[]{l}\xi_{L}=\frac{{2{l_{GPL}}}}{{{t_{GPL}}}},\;\xi_{W}=\frac{{2{w_{GPL}}}}{{{t_{GPL}}}},\;\eta_{L}=\frac{{\left({{E_{GPL}}/{E_{m}}}\right)-1}}{{\left({{E_{GPL}}/{E_{m}}}\right)+\xi_{L}}},\;\eta_{W}=\frac{{\left({{E_{GPL}}/{E_{m}}}\right)-1}}{{\left({{E_{GPL}}/{E_{m}}}\right)+\xi_{W}}},\end{array}$$

$${\rho }_1={\rho }_{GPL}{V}_{GPL}+{\rho }_m{V}_m,$$

$${\nu }_1={\nu }_{GPL}{V}_{GPL}+{\nu }_m{V}_m,$$

$$\mathbf{u}(x,y,z)={\mathbf{u}}^0(x,y)+z{\mathbf{u}}^1(x,y)+f(z){\mathbf{u}}^2(x,y),$$

$${\bf{u}}=\left\{{\begin{array}[]{*{20}{c}}u\\ v\\ w\end{array}}\right\},{\rm{}}{{\bf{u}}^{0}}=\left\{{\begin{array}[]{*{20}{c}}{{u_{0}}}\\ {{v_{0}}}\\ {{w_{0}}}\end{array}}\right\},\,{\rm{}}\,{{\bf{u}}^{1}}=-\left\{{\begin{array}[]{*{20}{c}}{{w_{0,x}}}\\ {{w_{0,y}}}\\ 0\end{array}}\right\},{\rm{}}{{\bf{u}}^{2}}=\left\{{\begin{array}[]{*{20}{c}}{{\theta_{x}}}\\ {{\theta_{y}}}\\ 0\end{array}}\right\},$$

$$w_{0,x}={\beta }_x,w_{0,y}={\beta }_y.$$

$${{\bf{u}}^{0}}={\left\{{\begin{array}[]{*{20}{c}}{{u_{0}}}\\ {{v_{0}}}\\ {{w_{0}}}\end{array}}\right\}},{{\bf{u}}^{1}}=-{\left\{{\begin{array}[]{*{20}{c}}{{\beta_{x}}}\\ {{\beta_{y}}}\\ 0\end{array}}\right\}},{{\bf{u}}^{2}}={\left\{{\begin{array}[]{*{20}{c}}{{\theta_{x}}}\\ {{\theta_{y}}}\\ 0\end{array}}\right\}.}$$

$${\epsilon }_{ij}=\frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}+\frac{\partial u_k}{\partial x_i}\frac{\partial u_k}{\partial x_j}\right).$$

$$\mathbf{\epsilon }={\left\{{\mathbf{\epsilon }}_{xx},{\mathbf{\epsilon }}_{yy},{\mathbf{\gamma }}_{xy}\right\}}^T={\mathbf{\epsilon }}_m+z{\mathbf{\kappa }}_1+f(z){\mathbf{\kappa }}_2,$$

$$\mathbf{\gamma }={\left\{{\gamma }_{xz},{\gamma }_{yz}\right\}}^T={\mathbf{\epsilon }}_s+f^{\prime }(z){\mathbf{\kappa }}_s,$$

$$\begin{array}[]{l}{{\boldsymbol{\varepsilon}}_{m}}=\left\{{\begin{array}[]{*{20}{c}}{{u_{0,x}}}\\ {{v_{0,y}}}\\ {{u_{0,y}}+{v_{0,x}}}\end{array}}\right\}+\frac{1}{2}\left\{{\begin{array}[]{*{20}{c}}{w_{0,x}^{2}}\\ {w_{0,y}^{2}}\\ {2w_{0,xy}}\end{array}}\right\}={\boldsymbol{\varepsilon}}_{m}^{L}+{\boldsymbol{\varepsilon}}_{m}^{NL},\\ {{\boldsymbol{\kappa}}_{1}}=-\left\{{\begin{array}[]{*{20}{c}}{{\beta_{x,x}}}\\ {{\beta_{y,y}}}\\ {{\beta_{x,y}}+{\beta_{y,x}}}\end{array}}\right\},{{\boldsymbol{\kappa}}_{2}}=\left\{{\begin{array}[]{*{20}{c}}{{\theta_{x,x}}}\\ {{\theta_{y,y}}}\\ {{\theta_{x,y}}+{\theta_{y,x}}}\end{array}}\right\},{\boldsymbol{\varepsilon}_{s}}=\left\{{\begin{array}[]{*{20}{c}}{{w_{0,x}}-{\beta_{x}}}\\ {{w_{0,y}}-{\beta_{y}}}\end{array}}\right\},{\boldsymbol{\kappa}_{s}}=\left\{{\begin{array}[]{*{20}{c}}{{\theta_{x}}}\\ {{\theta_{y}}}\end{array}}\right\}\end{array}$$

$${\boldsymbol{\varepsilon}}_{m}^{NL}=\frac{1}{2}\left[{\begin{array}[]{*{20}{c}}{{w_{0,x}}}&0\\ 0&{{w_{0,y}}}\\ {{w_{0,y}}}&{{w_{0,x}}}\end{array}}\right]\left\{{\begin{array}[]{*{20}{c}}{{w_{0,x}}}\\ {{w_{0,y}}}\end{array}}\right\}=\frac{1}{2}{\boldsymbol{\Theta}}{\bf{\Lambda}}.$$

$$\left\{{\begin{array}[]{*{20}{c}}{\boldsymbol{\sigma}}\\ {\bf{D}}\end{array}}\right\}=\left[{\begin{array}[]{*{20}{c}}{\bf{c}}&{-{{\bf{e}}^{T}}}\\ {\bf{e}}&{\bf{g}}\end{array}}\right]\left\{{\begin{array}[]{*{20}{c}}{\hat{\boldsymbol{\varepsilon}}}\\ {\bf{E}}\end{array}}\right\},$$

$$\mathbf{E}=-\mathrm{grad}\varphi =-\nabla \varphi .$$

$${\bf{c}}=\left[{\begin{array}[]{*{20}{c}}{\bf{A}}&{\bf{B}}&{\bf{N}}&{\bf{0}}&{\bf{0}}\\ {\bf{B}}&{\bf{M}}&{\bf{G}}&{\bf{0}}&{\bf{0}}\\ {\bf{N}}&{\bf{G}}&{\bf{H}}&{\bf{0}}&{\bf{0}}\\ {\bf{0}}&{\bf{0}}&{\bf{0}}&{{{\bf{A}}^{s}}}&{{{\bf{B}}^{s}}}\\ {\bf{0}}&{\bf{0}}&{\bf{0}}&{{{\bf{B}}^{s}}}&{{{\bf{D}}^{s}}}\end{array}}\right],$$

$$\begin{array}[]{l}\left({\bf A,\bf B,\bf M,\bf N,\bf G,\bf H}\right)=\int_{-h_{c}/2}^{h_{c}/2}{\left({1,z,{z^{2}},f(z),zf(z),{f^{2}}(z)}\right)}Q^{b}dz,\\ \left({{\bf A}^{s},{\bf B}^{s},{\bf D}^{s}}\right)=\int_{-h_{c}/2}^{h_{c}/2}{(1,f^{\prime}(z),{{f^{\prime}}^{2}}(z))}Q^{s}dz,\end{array}$$

$${\bf Q}^{b}=\frac{{{E(z)}}}{{1-\nu(z)^{2}}}\left[{\begin{array}[]{*{20}{c}}1&{{\nu(z)}}&0\\ {{\nu(z)}}&1&0\\ 0&0&{\frac{1-\nu(z)}{2}}\end{array}}\right],\;\;{\bf Q}^{s}=\frac{{{E(z)}}}{{2(1+\nu(z))}}\left[{\begin{array}[]{*{20}{c}}1&0\\ 0&1\end{array}}\right].$$

$${\bf{e}}=\left[{\begin{array}[]{*{20}{c}}0&0&0&0&{{e_{15}}}\\ 0&0&0&{{e_{15}}}&0\\ {{e_{31}}}&{{e_{32}}}&{{e_{33}}}&0&0\end{array}{\rm{}}\begin{array}[]{*{20}{c}}\end{array}}\right],{\bf{g}}=\left[{\begin{array}[]{*{20}{c}}{{p_{11}}}&0&0\\ 0&{{p_{22}}}&0\\ 0&0&{{p_{33}}}\end{array}}\right].$$

$$\delta {\int }_{t_1}^{t_2}Ldt=0,$$

$$L=\frac{1}{2}{\int }_{\Omega }\left(\rho {\dot{\mathbf{u}}}^T\dot{\mathbf{u}}-{\mathbf{\sigma }}^T\mathbf{\epsilon }+{\mathbf{D}}^T\mathbf{E}\right)\mathrm{d}\Omega +{\int }_{{\Gamma }_s}{\mathbf{u}}^T{\mathbf{f}}_s\mathrm{d}{\Gamma }_s-{\int }_{{\Gamma }_{\varphi }}\varphi {\mathbf{q}}_s\mathrm{d}{\Gamma }_{\varphi }+\sum {\mathbf{u}}^T{\mathbf{F}}_p-\sum \varphi {\mathbf{Q}}_p,$$

$$\begin{array}[]{l}\int_{{t_{1}}}^{{t_{2}}}{\int_{\Omega}{\left({\rho{\bf{\ddot{u}}}\delta{\bf{u}}^{T}-{\bf{\boldsymbol{\sigma}}}^{T}\delta{\boldsymbol{\varepsilon}}+{\bf{D}}^{T}\delta{\bf{E}}}\right)}}{\rm{d}}\Omega{\rm{dt+}}\int_{{t_{1}}}^{{t_{2}}}{\int_{{\Gamma_{s}}}{{{\bf{f}}_{s}}\delta{\bf{u}}^{T}}}{\rm{d}}{\Gamma_{s}}{\rm{dt-}}\int_{{t_{1}}}^{{t_{2}}}{\int_{{\Gamma_{\phi}}}{{{\bf{q}}_{s}}\delta{\phi}}}{\rm{d}}{\Gamma_{\phi}}{\rm{dt}}\\ \\ {\rm{+}}\int_{{t_{1}}}^{{t_{2}}}{\sum{{{\bf{F}}_{p}}}}\delta{{\bf{u}}^{T}}dt-\int_{{t_{1}}}^{{t_{2}}}{\sum{{{\bf{Q}}_{p}}}}\delta\phi dt=0.\end{array}$$

Full text

Highlights

Bézier extraction based IGA is successfully implemented for geometrically nonlinear static and dynamic analyses of smart FG plate structure.

• 2.

Two porosity distributions and three GPL dispersion patterns along the thickness direction is examined.

• 3.

Influence of porosity coefficient, weight fraction of GPLs and external electrical voltage are investigated.

• 4.

Active control of nonlinear dynamic responses via sensor and actuator layers is shown.

• 5.

The best reinforcement performance of piezoelectric FG porous plate structures is explored.

1 Introduction

A large demand for high structural performance and multifunctionality having superior mechanical properties and chemical stability in engineering applications has been significantly increased in the last few years. With cellular structures, porous materials having outstanding properties such as excellent energy absorption, lightweight, heat resistance has been broadly employed in various engineering fields [1, 2]. However, the existence of internal pores in metal matrix leads to a significant reduction in terms of structural stiffness [3]. To make up for this shortcoming, reinforcement with carbonaceous nanofillers e.g. carbon nanotubes (CNTs) [4, 5] and graphene nanoplatelets (GPLs) [6] into lightweight materials is considered as an excellent and practical choice to strengthen their mechanical properties. More importantly, this reinforcement also aims to preserve their potential for lightweight structures [7]. Compared to CNTs, GPLs have exhibited great potentials to become a good reinforcement candidate [8] since GPLs possess superior mechanical properties, a lower manufacturing cost, a larger specific surface area and two-dimensional geometry. To enhance the performance of structures, functionally graded (FG) porous structures reinforced by GPLs have been proposed in the literature to obtain the desired mechanical properties by adjusting the size and density of porosities as well as distribution of internal pores and GPLs [9, 10]. In terms of numerical analysis, numerous investigations have been carried out to examine the effect of internal pores and GPLs on structural behaviors under various conditions. Based on Timoshenko’s beam theory and Ritz method, the free vibration and buckling analyses of FG porous beams reinforced with GPLs are conducted by Kitipornchai et al. [11] while the nonlinear free vibration and postbuckling behaviors are reported by Chen et al. [12]. In addition, the nonlinear responses of arch and ring porous structures reinforced by GPLs are also examined [13, 14]. Yang et al. [15] employed the first-order shear deformation plate theory (FSDT) and Chebyshev-Ritz method to investigate the buckling and free vibration behaviors of FG porous plate structures with GPL reinforcement. Then, the nonlinear vibration and dynamic buckling of FG porous plates reinforced with GPLs resting on the elastic foundations are considered in [16, 17]. Shahverdi et al. [18] presented the postbuckling analysis for perfect/imperfect honeycomb core sandwich plate with GPL reinforcement. Based on isogeometric analysis (IGA), Li et al. [19] performed the static, free vibration and buckling analyses of FG porous plate structures reinforced by GPLs.

On the other hand, piezoelectric materials have been developed and widely employed in various engineering fields to produce smart structures. One of the essential and excellent features of this material type is the capacity transformation between mechanical and electrical energy, which is known as piezoelectric effect and converse phenomenon [20]. Regarding the analysis for plate structures embedded in piezoelectric layers, numerous researches have been carried out in order to predict their behaviors in the literature [21, 22, 23, 24]. Besides, FG carbon nanotubes reinforced composite plates (FG-CNTRC) embedded in piezoelectric layers also attracted remarkable attention of researchers. Alibeigloo [25, 26] presented the static and free vibration analyses of FG-CNTRC plate as well as cylindrical panel embedded in thin piezoelectric layers. By employing FSDT, Sharma et al. [27] reported the active vibration control of FG-CNTRC plates carrying piezoelectric sensor and actuator layers. Based on element-free IMLS-Ritz model and HSDT, Selim et al. [28, 29, 30] conducted a series of works regarding the active vibration control of FG plates reinforced by CNTs or GPLs embedded in piezoelectric layers. The dynamic response of laminated CNTRC plates incorporated with piezoelectric layers is addressed by Nguyen-Quang et al. [31] using IGA based on HSDT. In addition, Malekzade et al. [32] examined the free vibration of FG eccentric annular plates with GPLs reinforcement and embedded in piezoelectric layers. Most recently, the FG porous plate structures reinforced by GPLs and integrated piezoelectric layers are considered by Nguyen et al. [33] for static and dynamic analyses and Nguyen et al. [34] for active vibration control problems.

It is known that the different basic functions are applied for approximation of geometry models and solutions in the framework of traditional FEM which can lead to errors in computational process. In order to tackle this issue, IGA which employed non-uniform rational B-splines (NURBS) basis as shape functions are developed by Hughes et al. [35]. The main idea of IGA is fulfilled by using the same shape functions to describe the geometry model as well as to approximate the solution field. Till now, the IGA has been successfully applied to various engineering fields [36]. In comparison with standard FEM, the NURBS based IGA yields better accuracy and efficiency for many engineering issues, particularly for ones with complicated geometry models. The basic and review of IGA are presented in the established literature [37, 36]. Nevertheless, it should be emphasized that implementing NURBS based IGA approach seems not to easy as its basis functions are not confined to a unique element but instead span over the whole parametric domain. In addition, the calculation of stiffness matrix as well as force vector requires a mapping process between parametric and physical spaces. To overcome these obstacles, Borden et al. [38] first proposed Bézier extraction technique for NURBS by using Bernstein polynomials basis which is similar to the Lagrangian shape functions. Accordingly, Bézier extraction technique decomposes NURBS basis function into a linear combination of Bernstein polynomials basis and obtains $C^{0}$ continuity of isogeometric Bézier elements. By employing Bernstein polynomials like basis functions in Bézier extraction concept, the implementation of IGA becomes compatible with the traditional FEM while higher degree shape functions are efficiently exploited. As a result, the IGA approach can be easily embedded in most existing FEM codes while its advantages are still kept naturally. Several research works based on Bézier extraction concept can be found in the literature [39, 40, 41].

From the above literature review, it can be observed that the mentioned researches have only focused on plate structures integrated with piezoelectric layers in which core layers are often made of FGM or FG-CNTRC. Furthermore, the geometrically nonlinear analyses of piezoelectric FG plates under various dynamic loads are still somewhat limited. In the view of practical applications of high performance and smart engineering structures, the geometrically nonlinear static and dynamic responses of piezoelectric FG plate structures should be thoroughly examined to support the structural design process. In addition, the exploiting and integrating the outstanding features of IGA into the available FEM framework should be more considered in the literature because of the excellent computational effectiveness of them. This work aims to fill the research gap in the literature regarding the analysis and control of geometrically nonlinear responses of piezoelectric FG plates having the core layer composed of porous materials reinforced by GPLs. Accordingly, Bézier extraction technique based IGA in conjunction with $C^{0}$-HSDT is fully exploited in this study. By using the Bézier extraction, the IGA retains the element structure which allows the IGA approach to integrate conveniently into the existing FEM routine. The active control for geometrically nonlinear static and dynamic responses of FG plates with structural damping effect using a closed loop control via sensor and actuator layers is presented. In this study, the core layers are instituted by a combination of two porosity distributions and three GPL dispersion patterns across plate’s thickness, respectively while the piezoelectric layers are perfectly bonded on both top and bottom surfaces of core layer. The Newmark’s integration scheme in conjunction with Newton-Raphson technique is utilized for nonlinear analyses. Several numerical verifications are first conducted in order to demonstrate the accuracy and efficiency of the proposed approach. Numerically obtained results show that the responses of piezoelectric FG are significantly influenced by porosity coefficient, weight fraction of GPLs, distribution of porosity and GPLs into metal matrix as well as input voltages. More importantly, the best reinforcement performance for piezoelectric FG porous plate structures reinforced by GPLs is discovered in this study.

The outline of this work as follows. Section 2 provides theoretical formulations including material models, $C^{0}$-type HSDT, constitutive equations and weak form. The brief on Bézier extraction based IGA and discrete system equations are given in Section 3. The active control algorithm is described in Section 4. Section 5 presents numerical investigations for the geometrically nonlinear static and dynamic analyses along with active control. The paper is closed with some concluding remarks which can be found in Section 6.

2 Theoretical formulations

2.1 Material models

Consider a FG plate model as a sandwich structure whose core layer is porous material reinforced by GPLs and piezoelectric layers considered as face layers, as depicted in Fig. 1. The geometric dimensions of FG plate are taken as length $a$, width $b$ and total thickness $h=h_{c}+2h_{p}$, where $h_{c}$ and $h_{p}$ represents the thicknesses of core layer and each piezoelectric layer, respectively. In this work, the core layers are instituted by a combination, respectively, of two porosity distributions (symmetrical and asymmetrical) and three GPL dispersion patterns (symmetrical, asymmetrical and uniform) through the thickness direction, as shown in Fig. 3. The effective Young’s modulus $E(z)$, shear modulus $G(z)$ and mass density $\rho(z)$ for two porosity distributions can be mathematically expressed as follows

[TABLE]

in which porosity distribution functions $\lambda(z)$ are defined as

[TABLE]

where $E_{1}$ and $\rho_{1}$ are the maximum values of Young’s modulus and mass density in the thickness direction of core layer, respectively. The coefficient of porosity $e_{0}$ in Eq. (1) is given by

[TABLE]

in which $E^{\prime}_{1}$ and $E^{\prime}_{2}$ denote, respectively, the maximum and minimum values of Young’s modulus for core layers without GPLs reinforcement, as shown in Fig. 3(a). By employing the Gaussian Random Field (GRF) scheme [42], the mechanical properties of closed-cell cellular solids are determined as

[TABLE]

From Eq. (4), the mass density parameter $e_{m}$ in Eq. (1) are determined as

[TABLE]

In addition, the effective Poisson’s ratio $\nu(z)$ can be given as [43]

[TABLE]

where $\nu_{1}$ denotes the Poisson’s ratio of metal without internal pores. The parameter $\varsigma$ in Eq. (6) can be derived as follows

[TABLE]

Regarding the distribution of GPLs into metal matrix, the varied volume fraction along the $z$ direction for three GPL dispersion patterns, as illustrated in Fig. 3(b), can be formulated as follows

[TABLE]

where $S_{i1}$,$S_{i2}$ and $S_{i3}$ are the peak values of GPL volume fraction, in which $i=1,2$ corresponding to two different porosity distributions. Additionally, the relation of volume fraction and weight fraction of GPLs can be expressed as [33]

[TABLE]

From Eq. (9), the peak values of GPL volume fraction in Eq. (8) can be defined. Besides, based on Halpin-Tsai micromechanics model [44], the value of Young’s modulus $E_{1}$ of core layer with GPLs reinforcement without internal pores is defined as

[TABLE]

in which

[TABLE]

in which $w_{GPL}$, $l_{GPL}$ and $t_{GPL}$ are fundamental dimensions of GPLs including the average width, length and thickness, respectively. It should be mentioned that the mechanical properties of metal matrix and GPLs are symbolized with subscript symbols $m$ and $GPL$, respectively. Finally, the mass density $\rho_{1}$ and Poisson’s ratio $\nu_{1}$ can be determined based on the rule of mixture [45]

[TABLE]

[TABLE]

in which the relation of $V_{GPL}$ and $V_{m}$ can be given as $V_{GPL}+V_{m}=1$.

2.2 $C^{0}$-type higher-order shear deformation theory

Considering a plate carrying a domain ${\bf V}=\Omega\times(\frac{-h}{2},\frac{h}{2})$, in which $\Omega\in\mathbb{R}^{2}$. Assuming that the displacement components of plate in the $x$, $y$, $z$ directions and rotations about the $y$- and the $x$- axes are $u_{0},v_{0},w_{0},\theta_{x}$ and $\theta_{y}$, respectively. Based on the higher-order shear deformation theory [46], the displacement field at an any point in plate can be described as follows

[TABLE]

where

[TABLE]

where the subscript symbols $x$ and $y$ denote the derivative of any function corresponding to $x$ and $y$ directions, respectively. Meanwhile, the function $f(z)$ describes the shear strain and stress distribution along the thickness, as found in [47]. In this work, the famous third-order function proposed by Reddy is utilized as $f(z)=z-\frac{4z^{3}}{3h^{2}}$ [48].

In order to avoid the high-order derivations in approximate formulation and enforce boundary conditions conveniently, two extra assumptions are formulated as

[TABLE]

By substituting Eq. (16) into Eq. (15), yields

[TABLE]

It can be seen that the strain fields in Eq. (17) only involve the $C^{0}$-continuity of approximate field. The Green strain vector of a bending plate can be written in compact form as follows

[TABLE]

By employing the von Kármán assumptions, the relationship of strain and displacement is given as

[TABLE]

where

[TABLE]

where the nonlinear strain component is expressed as follows

[TABLE]

2.3 Constitutive equations

In this study, the constitutive relationship of FG plate embedded in piezoelectric layers is expressed as [49]

[TABLE]

in which ${\boldsymbol{\sigma}}$ and $\hat{\boldsymbol{\varepsilon}}={{\rm{\{}}{{\boldsymbol{\varepsilon}}},\;{{\boldsymbol{\gamma}}}{\rm{\}}}^{T}}$ represent, respectively, the stress and strain vectors of mechanical field. $\bf D$ stands for dielectric displacement vector; Meanwhile, $\bf c$, $\bf e$ and $\bf g$ are the mechanical, piezoelectric and dielectric constant matrices, respectively. In addition, the electric field vector $\bf E$ is calculated based on the electric potential field $\phi$ as follows [50]

[TABLE]

in which $\phi$ stands for the electric potential difference through the piezoelectric layer. The material constant matrix $\bf c$ can be expressed as

[TABLE]

in which

[TABLE]

in which

[TABLE]

where $E(z)$ and $\nu_{(}z)$ are the effective Young’s modulus and Poisson’s ratio, respectively. Additionally, the constant matrices of piezoelectric material, as presented in Eq. (22), including the stress piezoelectric $\bf e$ and the dielectric $\bf g$ are given by [51]

[TABLE]

2.4 Weak form of governing equations

The governing equations of motion for piezoelectric FG plate can obtain by applying the Hamilton’s variational principle which are given as follows [52]

[TABLE]

in which $t_{1}$ and $t_{2}$ represent the starting and end time values, respectively. The general energy function $L$ is the summation of kinetic energy, strain energy, dielectric energy and external work as follows

[TABLE]

where $\rho$ represents the mass density; $\bf u$ and $\dot{\bf u}$ denote the mechanical displacement and velocity vectors, respectively; ${\bf f}_{s}$ and ${\bf F}_{p}$ are, respectively, the external mechanical surface and concentrated load vectors while ${\bf q}_{s}$ and ${\bf Q}_{p}$ are the external surface and point charges, respectively; $\Gamma_{s}$ and $\Gamma_{\phi}$ indicate the external mechanical and the electrical loading surface, respectively.

Then, the variational form of Eq. (29) can be expressed as follows

[TABLE]

3 Isogeometric analysis based on Bézier extraction

3.1 B-spline and NURBS basis functions

In the parametric space for one dimensional (1D) problem, B-spline basis functions can be established by non-decreasing set real-valued coordinates called knot vector ${\bf{\Xi}}=\left\{{{\xi_{1}},{\xi_{2}},...,{\xi_{n+p+1}}}\right\}$, in which $\xi_{i}$ is the $i$th knot with $i=1,...n+p+1$ indicates the knot index. Meanwhile, $n$ and $p$ are the number of univariate B-spline basis functions and its order, respectively. For a given knot vector $\bf{\Xi}$, the B-spline basis functions are defined according to recursive form [35, 37]

[TABLE]

[TABLE]

Then, B-spline curves can be constructed by taking a linear combination of B-spline basis functions and control points ${P_{i}}\left({i=1,2,...,n}\right)$ as

[TABLE]

The B-spline basis functions in two dimensional (2D) space can also obtain by taking a tensor product of two 1D B-splines basis functions. The B-splines surfaces are then expressed as follows [35]

[TABLE]

in which ${N_{i,p}}$ and ${M_{j,q}}$ represent the 1D B-splines basis functions having orders $p$ and $q$ in the $\xi$ and $\eta$ directions corresponding with two knot vectors ${\bf{\Xi}}=\left\{{{\xi_{1}},{\xi_{2}},...,{\xi_{n+p+1}}}\right\}$ and ${\boldsymbol{\Theta}}=\left\{{{\eta_{1}},{\eta_{2}},...,{\eta_{m+q+1}}}\right\}$, respectively.

In order to describe exactly some conic shapes such as circle, cylinder, sphere and ellipsoid, the NURBS is introduced based on the B-spline and a set of weights [37]. Accordingly, the NURBS basis functions can be expressed as

[TABLE]

where ${w_{i,j}}$ represents the weight values. Finally, the NURBS surfaces can be determined as follows

[TABLE]

3.2 Bézier extraction of NURBS

The main aim of Bézier extraction technique is to instead the NURBS basis functions by the $C^{0}$-continuous Bernstein polynomial basis functions described over Bézier elements which have the similar element structure with standard FEM. As a result, the IGA is straightforwardly performed as well as can be integrated into most available FEM framework. It is well known that the B-spline basis functions of $p$th order possess $C^{p-k}$ continuity across each interior knots, in which $k$ represents the multiplicity of knots in the knot value. Hence, we can obtain the $C^{0}$-continuity by inserting the new knots into the original knot vector until $k$ is equal to $p$. By inserting a new knot $\bar{\xi}\in\left[{{\xi_{k}},{\xi_{k+1}}}\right]$ with $(k>p)$ into the original knot vector ${\bf{\Xi}}=\left\{{{\xi_{1}},{\xi_{2}},...,{\xi_{n+p+1}}}\right\}$, a new set of control points ${\bar{\bf P}_{i}}$ are found and expressed as follows [35]

[TABLE]

in which

[TABLE]

Then, the Bézier extraction operator can be established by using the new set of knots as follows [38, 53]

[TABLE]

By employing the Bézier extraction operator ${\bf C}^{j}$, a new Bézier control points ${\bf P}^{b}$ associated with Bernstein polynomial basis can obtain as follows [54]

[TABLE]

in which the whole Bézier extraction operator $\bf C$ can be given as

[TABLE]

It should be noted that the above insertion into the original knot vector does not lead to any change of the geometry. Consequently, the B-spline surface can also obtain by using Bézier control points and Bernstein polynomials as

[TABLE]

where 2D Bernstein polynomials ${\bf B}\left({\xi,\eta}\right)$ are defined recursively as

[TABLE]

in which

[TABLE]

Comparing Eqs. (34) and (42), yields

[TABLE]

By substituting Eq. (40) into Eq. (45), the relationship between the B-spline and Bernstein polynomial functions can obtain as follows

[TABLE]

Based on Eq. (46), the NURBS basis functions can be expressed by Bernstein polynomials as follows

[TABLE]

where $\bf{W}$ denotes the local NURBS weights in diagonal matrix. Meanwhile, the weight functions $W\left({\xi,\eta}\right)$ are described based on Bernstein basis functions as follows

[TABLE]

in which $\bf{w}$ and $\bf{w}^{b}$ indicate, respectively the weights of the NURBS and Bézier. In addition, the relationship of Bézier and NURBS control points can be given as

[TABLE]

3.3 Discrete system equations

3.3.1 Mechanical displacement field

By employing Bézier extraction of NURBS, the approximation of mechanical displacement field $\bf{u{(\xi,\eta)}}$ for FG plate is expressed as

[TABLE]

in which $m\times n$ represents the number of basis functions and ${R}_{I}^{e}\left({\xi,\eta}\right)$ denotes the NURBS basis functions which are presented in Eq. (47). The vector {\bf d}_{I}={\left\{{\begin{array}[]{*{20}{c}}{{u_{0I}}}&{{v_{0I}}}&{{w_{0I}}}&{{\beta_{xI}}}&{{\beta_{yI}}}&{{\theta_{xI}}}&{{\theta_{yI}}}\end{array}}\right\}^{T}} includes the nodal degrees of freedom associated with control point I.

Substituting Eq. (50) into Eq. (20), the strain components are expressed as

[TABLE]

in which {\bf{B}}_{I}^{L}={\left[{\begin{array}[]{*{20}{c}}{{\bf{B}}_{I}^{m}}&{{\bf{B}}_{I}^{b1}}&{{\bf{B}}_{I}^{b2}}&{{\bf{B}}_{I}^{s1}}&{{\bf{B}}_{I}^{s2}}\end{array}}\right]^{T}} is given as follows

[TABLE]

and ${\bf{B}}_{I}^{NL}={\boldsymbol{\Theta}}{\boldsymbol{\Lambda}}$ with

[TABLE]

3.3.2 Electric potential field

The electric potential field on piezoelectric layers is estimated by dividing each layer into $n_{sub}$ sublayers through the thickness. In each sublayer, the electric potential variation is assumed to be linear across thickness and calculated as follows [51]

[TABLE]

where ${\bf R}_{\phi}$ indicates the shape function of the electric potential function which is determined in Eq. (47) with order $p=1$. The function {\boldsymbol{\phi}}_{I}=\left\{{\begin{array}[]{*{20}{c}}{{\phi^{i-1}}},&{{\phi^{i}}}\end{array}}\right\} with $i=1,2,...,{n_{sub}}$ is the electric potentials at the upper and lower surfaces of each sublayer.

Assuming that at the same height, the values of electric potentials are set to be equal [49]. As a result, the electric potential field $\bf E$ for each sublayer element can be given as

[TABLE]

where

[TABLE]

3.4 Governing equation of motion

By substituting Eqs. (50) and (55) into Eq. (30), the final form of governing equation of motion can obtain as follows [51]

[TABLE]

where

[TABLE]

in which

[TABLE]

where

[TABLE]

and

[TABLE]

where ${I}_{i}$ with $i=1:6$ are defined as

[TABLE]

It should be noted that the electric potential field $\bf E$ depends only on following to $z$ direction leads to the matrix ${\bf K}_{u\phi}$ in Eq. (58) can be rewritten as

[TABLE]

Finally, by substituting the second equation into the first one of Eq. (57), yields

[TABLE]

4 Active control analysis

Consider a piezoelectric FG plate, as illustrated in Fig. 4, for the active control static and dynamic responses. While the lower layer is a piezoelectric sensor named with the subscript $s$, the upper one represents a piezoelectric actuator labeled with the subscript $a$. Because of the piezoelectric effect, the electric charge can be generated and collected in the sensor layer when the structures are deformed. After that, through an appropriate electronic circuit, this electric charge is amplified and converted into the voltage signal before being sent and applied to the actuator layer via closed loop control algorithm. Due to the converse piezoelectric phenomenon, the corresponding deformations can occur in structures which lead to damping the structural dynamic responses.

By ignoring the external charge $\bf Q$, the generated potential electric in sensor layer can be estimated based on the second equation of Eq. (57)

[TABLE]

Then, the sensor charge resulted because of the corresponding deformation can be expressed as

[TABLE]

In this study, two constant control gains of the displacement $G_{d}$ and velocity feedback $G_{v}$, respectively, are employed to couple the input and output voltage vectors as follows [55]

[TABLE]

By substituting Eqs. (67) and (65) into the second equation in Eq. (57), the magnitude of actuator layer charge can be archived as follows

[TABLE]

Next, by substituting Eq. (68) into Eq. (64), one obtains

[TABLE]

in which

[TABLE]

and the active damping matrix $\bf C$ is expressed as

[TABLE]

Considering the structural damping effect, Eq. (69) can be rewritten as

[TABLE]

where the Rayleigh damping matrix ${\bf C}_{R}$ is determined based on a linear combination between $\bf M$ and ${\bf K}_{uu}$, as follows

[TABLE]

in which $\alpha_{R}$ and $\beta_{R}$ denote Rayleigh damping coefficients which can obtain from experiments.

5 Numerical examples

In this work, Newton-Rapshon iterative procedure is employed to obtain the solution of nonlinear problems. For the geometrically nonlinear dynamic analysis of FG plate under various dynamic loads, which the equations depend on both the time domain and variable displacement, Newmark’s integration scheme is exploited. In all numerical examples, the PZT-G1195N piezoelectric is employed as well as perfectly glued on the upper and lower surfaces of plate and neglected the adhesive layers. In addition, a mesh of $15\times 15$ quadratic Bézier elements, as illustrated in Fig. 5, is utilized in order to model all square plate structures in this work.

5.1 Validation analysis

In this section, various numerical studies regarding the geometrically nonlinear static and dynamic analyses of isotropic as well as piezoelectric FG plates are carried out in order to perform the accuracy of the proposed approach. For geometrically nonlinear analysis, an isotropic square plate with a fully clamped (CCCC) boundary conditions is considered to show the validity of the present formulation. The plate structure is subjected to a uniform transverse load while the width-to-thickness ratio $(a/h)$ is taken equal to $100$. The material properties are $E=3\times 10^{7}$ psi and $\nu=0.316$. In this example, the normalized central deflection and load parameter are given as $\overline{w}=w/h$ and $P=q_{0}a^{4}/(Eh^{2})$, respectively. Table LABEL:tab:ISO_CCCC presents the normalized central deflections which are compared with those of Levy’s analytical solution [56], Urthaler and Reddy’s mixed FEM using FSDT [57] and Nguyen et al. based on IGA and refined plate theory [58]. As can be observed that the present results are in good agreement with the existing analytical solution as well as other approximate results.

Next, in order to verify the accuracy of the proposed approach for the geometrically nonlinear dynamic analysis, a fully simply supported (SSSS) boundary conditions orthotropic square plate subjected to a uniform load with $q_{0}=1.0$ MPa is conducted. The material properties and the geometry dimensions are given as follows: Young’s modulus $E_{1}=525$ GPa, $E_{2}=21$ GPa, shear modulus $G_{12}=G_{23}=G_{13}=10.5$ GPa, Poisson’s ratio $\nu=0.25$, mass density $\rho=800$ kg/m3, length $a=b=250$ mm and thickness $h=5$ mm. Fig. 6 depicts the geometrically nonlinear dynamic response of square plate. It can be observed that the present result matchs very well with that obtained from the finite strip method, as reported by Chen at al. [59].

Regarding the plate structures integrated with piezoelectric layers, a cantilever piezoelectric FG square plate is thoroughly considered to verify the accuracy of the present approach. The FG plate which is bonded by two piezoelectric layers on both the lower and upper surfaces is made of aluminum oxide and Ti-6Al-4V materials and is subjected to simultaneously a uniform transverse load with $q_{0}$=100 N/m2 and various input voltage values. The mechanical properties of two materials are given in Table 2. In this study, the rule of mixture [45] is employed to calculate the effective material properties of FG core layer. The FG plate has a side length $a=b=0.4$ m whilst the thickness of core layer and each piezoelectric layer are $h_{c}=5$ mm and $h_{p}=0.1$ mm, respectively. The centerline linear deflections of piezoelectric FG square plate are plotted in Fig. 7 while the tip node deflections are also listed in Table LABEL:tab:Deflec_Tip with various material index $n$. The obtained results are compared with those reported in [60] based on a cell-based smoothed discrete shear gap method (CS-DSG3) using FSDT. It can be observed that the present results generally agree well with the reference solutions.

In the next Sections, investigations into the geometrically nonlinear static and dynamic responses of piezoelectric FG porous plate with GPLs reinforcement will be presented.

5.2 Geometrically nonlinear static analysis

The geometrically nonlinear static analysis of a piezoelectric FG plate subjected to uniform load with $\overline{q}=q_{0}\times 10^{3}$ is first addressed. A SSSS piezoelectric FG square plate made of aluminum oxide and Ti-6Al-4V has $a=b=0.2$ m, $h_{c}=2.0$ mm and $h_{p}=0.1$ mm. Fig. 8 illustrates the influence of material index $n$ on the normalized linear and nonlinear central deflections of piezoelectric FG plates under mechanical load. As can be observed that with an increase of material index $n$ leads to a decrease gradually of the deflection. The largest deflection occurs with $n=0.0$ where the plate includes only of Ti-6Al-4V leads to a decrease in the bending stiffness. Furthermore, the values of the central deflection for geometrically nonlinear analysis are always smaller than that of the linear one. This difference would gradually reduce with the increase of material index.

In the next example, we consider a SSSS piezoelectric FG plate having core layer constituted by a combination of two porosity distribution types and three GPL patterns, respectively. The FG plate is subjected to sinusoidally distributed load defined as $q=q_{0}sin(\pi x/a)sin(\pi y/b)$ in which $q_{0}=1.0$ MPa. The geometry dimensions of plate are taken as $a=b=0.4$ m, $h_{c}=20$ mm and $h_{p}=1.0$ mm. In this study, the copper is chosen as the metal matrix whose material properties are given in Table 2. Meanwhile, the basic dimensions of GPLs are ${l_{GPL}}=2.5\;\mu m$, ${w_{GPL}}=1.5\;\mu m$, ${t_{GPL}}=1.5\;nm$. Fig. 9 examines the influence of porosity coefficient on the nonlinear normalized deflection of piezoelectric FG porous plate with pattern $A$ $(\Lambda_{GPL}=1.0\;wt.\;\%)$ for two porosity distributions, respectively. It can be observed that an increase of the porosity coefficient leads to the corresponding increase of the nonlinear deflection since the higher density of internal pores in material yields the reduction in the stiffness of structures.

In addition, Fig. 10 depicts the influence of weight fraction and GPL dispersion patterns on the behaviors of FG porous plate with $e_{0}=0.2$ and two porosity distributions, respectively. We observe that the effective stiffness of FG porous core layer is significantly strengthened after adding a small amount of GPLs ($\Lambda_{GPL}=1.0\;wt.\;\%)$ into metal matrix as evidenced by a decrease of the nonlinear normalized deflection. More importantly, the reinforcement effect of GPLs depends significantly on the dispersion type of GPLs into metal matrix. Accordingly, with the same weight fraction of GPLs, the pattern $A$, where GPLs are dispersed symmetrically through the midplane of porous core layer, achieves the smallest nonlinear deflection while the asymmetric pattern $B$ provides the largest one. For further illustration, Fig. 11 depicts the variation of nonlinear normalized deflection of piezoelectric FG porous plate with GPLs reinforcement which is constituted by a combination of porosity distribution 1 and three different GPL dispersion patterns, respectively with load parameter 10.

The combined influences of two porosity distributions and three GPL patterns on the nonlinear deflection of FG plate structure with $\Lambda_{GPL}=1.0\;wt.\;\%$ and $e_{0}=0.4$ is thoroughly examined. As evidently depicted in Fig. 12, for all the considered cases, the combination of porosity distribution 1 and pattern $A$ always obtains the best reinforcement performance in the geometrically nonlinear static problems. This indicated that the FG plate structures where the internal pores are dispersed on the midplane and GPLs are distributed around the bottom and top surfaces yield the optimum reinforcement.

5.3 Geometrically nonlinear dynamic analysis

In this part, the geometrically nonlinear dynamic behaviors of a CCCC piezoelectric FG porous plate with GPL reinforcement are investigated. The basic dimensions as well as the material properties of FG plate are the same previous example. The FG plate is assumed to be subjected to time-dependent sinusoidally distributed transverse loads expressed as follows $q=q_{0}sin(\pi x/a)sin(\pi y/b)F(t)$. In this work, four different dynamic loads are examined via the function $F(t)$ as follows

[TABLE]

in which $q_{0}=100$ MPa, $\gamma=330s^{-1}$. Fig. 13 illustrates the time history $F(t)$ with $t_{1}=3\times 10^{-3}$ s.

Fig. 14 illustrates the influence of porosity coefficient on the nonlinear dynamic response of piezoelectric FG porous plate with porosity distribution 1 and pattern $A$ $(\Lambda_{GPL}=1.0\;wt.\;\%)$ under step and sinusoidal loads. It can be observed that when the porosity coefficient increases, the amplitude of the transverse deflection increases as well while the period of motion does not seem to affect. This is because the presence of porosities in core layer leads to reduction the capacity of itself against external excitation. Furthermore, Fig. 15 performs the influence of weight fraction and GPL dispersion pattern on the nonlinear dynamic response of piezoelectric FG porous plate with $e_{0}=0.2$ and porosity distribution 2 corresponding to triangular and explosive blast loads. As expected, the smaller magnitude of the deflection can obtain when the weight fraction of GPLs into metal matrix increases. Again, the dispersion of GPLs into metal matrix significantly influence to the reinforcement performance of structure that dispersion pattern $A$ provides the smallest magnitude of the deflection.

Next, the combination influences of various porosity distribution types and GPL dispersion patterns on the nonlinear dynamic behaviour is also examined and indicated in Fig. 16. For this specific example, the porous core layer has porosity coefficient $e_{0}=0.4$ and GPL weight fraction $\Lambda_{GPL}=1.0\;wt.\;\%$. As clearly demonstrated in Fig. 16, the combination of porosity distribution 1 and GPL dispersion pattern $A$ always provides the best reinforcement as evidenced by obtaining the smallest amplitude of the deflection. Moreover, the linear and nonlinear dynamic responses of FG porous plate with porosity distribution 2 $(e_{0}=0.3)$ and GPL dispersion pattern $C$ ($\Lambda_{GPL}=1.0\;wt.\;\%$) under triangular and sinusoidal loads are also plotted in Fig. 17. As can be observed, the geometrically nonlinear responses generally obtain smaller magnitudes of the deflection and periods of motion.

5.4 Static and dynamic responses active control

In this section, the active control for static and dynamic behaviors of FG porous plate with GPL reinforcement via sensor and actuator layers is performed. First, to demonstrate the accuracy of the present approach, the active control for linear static behaviors of a SSSS FG plate subjected to a uniform transverse load with $q_{0}$=100 N/m2 is examined. The FG plate having $a=b=0.2$ m, $h_{c}=1$ mm and $h_{p}=0.1$ mm are composed of Ti-6Al-4V and aluminum oxide materials with material index $n=2.0$. Fig. 18 illustrates the linear static deflections of FG plate with various displacement feedback control gains $G_{d}$. As can be observed that the obtained results from the proposed approach coincide well with the reference solution [60] employed the CS-DSG3 based on FSDT. As expected, when the displacement feedback control gain $G_{d}$ increases, the linear static deflection of FG plate decreases. Furthermore, by using a constant gain velocity feedback $G_{v}$ and a closed loop control, the active control for the linear dynamic responses of FG plate is also carried out in this work. Accordingly, we consider a FG plate initially subjected to a uniform transverse load $q_{0}$=100 N/m2 and then suddenly eliminated. It should be noted that to save the computational cost, the modal superposition is exploited and the first six modes are only employed in the modal space analysis. Additionally, for each mode the initial modal damping ratio is assumed to be 0.8 $\%$. Fig. 19 illustrates the linear dynamic behavior of central deflection for FG plate corresponding to two values of $G_{v}$. The results which are generated from present method agree well with the reference solution [60].

Next, the active control for nonlinear static responses of a SSSS FG porous plate reinforced with GPLs is further investigated. The core layer of FG plate consisting of combined of porosity distribution 1 and pattern A, which provides the best structural performance, is considered in this example. The material properties are the same in Section 5.2 while the geometric dimensions are given as $a=b=0.4$ m, $h_{c}=20$ mm and $h_{p}=1$ mm. The FG plate is subjected to a sinusoidally distributed load expressed as $q=q_{0}sin(\pi x/a)sin(\pi y/b)$ with $q_{0}=1.0$ MPa. Fig. 20 depicts the nonlinear normalized central deflection with $e_{0}=0.4$ and $\Lambda_{GPL}=1.0\;wt.\;\%$ corresponding to various displacement feedback control gains $G_{d}$. As can be seen that the normalized central deflection of FG porous plate decreases significantly when the values of displacement feedback control gain increase.

The last example aims to show the active control for geometrically nonlinear dynamic behaviors of a CCCC FG porous plate with GPL reinforcement under sinusoidally distributed transverse loads. The square plate with the porosity distribution 1 $(e_{0}=0.4)$ and dispersion pattern $A$ $(\Lambda_{GPL}=1.0\;wt.\;\%)$ has $a=b=0.2$ m, $h_{c}=10$ mm and $h_{p}=0.1$ mm. Fig. 21 illustrates the nonlinear dynamic responses of normalized central deflection corresponding to various dynamic load types as well as the velocity feedback control gains $G_{v}$. We observe that in case of uncontrolled ($G_{v}=0$), the nonlinear dynamic behaviors of FG porous plate still attenuate with respect to time since the structural damping effect is taken into account in this study. More importantly, the amplitude of central deflection can be suppressed more faster in the case controlled via higher velocity feedback control gain values. As a result, depending on a specific case, the nonlinear dynamic responses of structure such as deflection, oscillation time or even both can be adjusted in order to satisfy an expectation through designing an appropriate value for $G_{v}$. In addition, Fig. 22 shows the influence of velocity feedback control gain $G_{v}$ on the linear and nonlinear dynamic responses of plate subjected to step load. As would be expected, the nonlinear dynamic responses obtain smaller magnitudes of central deflection and periods of motion. However, it should be emphasized that the values of feedback control gain could not be enlarged without limit as piezoelectric materials have their own failure voltage values.

6 Conclusions

In this study, the IGA based on Bézier extraction of NURBS and $C^{0}$-HSDT was successfully implemented for the geometrically nonlinear static and dynamic analyses of piezoelectric FG porous plates with GPLs reinforcement. The motion equations were derived based on the $C^{0}$-HSDT in conjunction with von Kármám strain assumptions. While the mechanical displacements were estimated by using $C^{0}$-HSDT based on Bézier extraction of NURBS, the electric potential field was assumed as a linear function through the thickness of each piezoelectric layer. Two porosity distributions and three dispersion patterns of GPLs with various related parameters were thoroughly examined via several numerical investigations. The control algorithms based on the constant displacement and velocity feedback were utilized to control the nonlinear static and dynamic responses of FG porous plate structures reinforced with GPLs. A number of key remarks can be drawn from the present work:

By applying Bernstein polynomials as basis functions in Bézier extraction, the IGA can easily be integrated into most existing FEM structures while its advantages are still maintained effectively.

• 2.

The stiffness of plate structures is significantly enhanced after adding a small amount of GPLs into metal matrix, while the existence of internal pores yields the reduction in both structural weight and reinforcement effect.

• 3.

The distribution of porosities and GPLs in metal matrix plays a crucial role in reinforcing performance. A combination of porosity distribution 1 with internal pores distributed around the midplane and GPL pattern $A$, where GPLs are close dispersed the top and bottom surfaces, achieve the best reinforcement performance.

• 4.

The nonlinear dynamic responses can be actively controlled through the velocity feedback control algorithm based on a closed-loop control.

Acknowledgements

This research was supported by the Grant (NRF-2017R1A4A1015660) from NRF (National Research Foundation of Korea) funded by MEST (Ministry of Education and Science Technology) of Korean government. The third author would like to thank the support of the RISE-project BESTOFRAC (734370)–H2020.