Machine learning assisted coupled frequency analysis of skewed multi-phase magnetoelectric composite plates with temperature and moisture dependent properties

Abstract

In this article, the application of an artificial neural network (ANN)-based machine learning (ML) strategy to predict the coupled frequency of geometrically skewed multiphase magnetoelectric (MME) composite plate exposed to hygrothermal environment is presented. The ANN model is trained using a dataset comprising more than one million simulations conducted using an in-house developed finite element formulation. The underlying multiphysics governing equations are derived using Hamilton’s principle and higher-order shear deformation theory (HSDT). The influence of the hygrothermal environment on the elastic stiffness of MME composites is defined by the empirical constants in the constitutive relations. Four prominent combinations of the empirical constants leading to different elastic stiffness relations have been considered in this study. Alongside, the influence of geometrical skewness on the coupled fundamental frequency is also assessed. For the training of the ANN model, the Levenberg–Marquardt optimization algorithm with 30 neurons along with a tangent sigmoid activation function is used. The trained ANN model is tested on an unseen dataset, different from the training data, and it is shown to accurately predict the natural frequency of MME plate with a maximum error of 2.3%. Further, excluding the training time and considering the computational time alone, the ANN model is found to be 6.3 times faster than the FE simulations. It is anticipated that such ML-based reduced order models can be effective in the design process, especially in complex multiphysics problems, such as the one considered in the work, involving a multitude of geometric, loading and material parameters.

Keywords:

• Artificial neural network
• machine learning
• magnetoelectric
• finite element
• hygrothermal
• empirical constant
• skew

Nomenclature
u, v, and w	=	displacement components in x-, y- and z- direction
u0, v0, and w0	=	mid-plane displacement components
a, b, h	=	ME plate’s length, width and thickness
Vf	=	volume fraction
Tp, Tk	=	potential and kinetic energy
θx and θy	=	rotations of normal in xz-plane and yz-plane
νx, νy	=	higher order rotational terms
φ	=	skew angle
εb, εs	=	bending and shear strains
γ, β	=	empirical constants
${\mathrm{\Omega } }^{\mathit{plate}}$	=	volume of the plate
ρ	=	density
ΔT, Δm	=	temperature and moisture gradient
[C0], [C]	=	temperature and moisture independent and dependent elastic stiffness coefficients, respectively
[Btb], [Brb], [Bts], [Brs]	=	strain-displacement matrices
[e], [q], [m], [η] and [µ]	=	piezoelectric, magnetostrictive, magnetoelectric, dielectric and permeability matrices
{σ}, {D} and {B}	=	stress, electric displacement, and magnetic flux density vectors
{E}, {H}	=	electric and magnetic field intensity
[Tt], [Tr]	=	transformation matrices for translational and rotational DoFs
[T1], [T2]	=	transformation matrices for stiffness matrices

1. Introduction

Composite structures are being extensively used for various engineering applications ranging from civil to aerospace structures, owing to their flexibility and tailorability. While their exceptional load-bearing capability is well-utilised by industry, recent research efforts are being focussed on imparting multifunctionalities into such composite structures which involves adding piezoelectric, piezomagnetic and electrochemical characteristics [1]. Multifunctional composites have shown promise in a range of applications such as energy harvesting [2–8], active vibration control [9–12], morphing [13–15] and energy storage [16], to name a few. Amongst different forms of multifunctional composites, magnetoelectric (ME) composites display a beneficial interaction among different fields (electric, elastic, and magnetic) which makes them a potential candidate for various smart structural applications such as energy harvesters [17–19] and vibration control [20, 21]. Alongside, the coupling exhibited by ME composites are predominantly high as opposed to other categories of smart multifunctional materials. The degree of coupling can be tweaked by varying the weight percentage of piezoelectric (PE) and piezomagnetic (PM) phases. The ME composites with different volume fractions of PE and PM phases are known as multiphase ME (MME) composites.

Multiple research studies exist to quantify the multiphysics response of ME composite structures under mechanical and other loading conditions. Unlike conventional composites where the mechanical properties are of primary interest, analyzing multifunctional composite structures can be increasingly complex due to their coupled behavior [22–24]. Further, exposure to hazardous environmental conditions can significantly influence the material properties and the structural responses, an example being the effect of hygrothermal (HT) conditions. Numerous researchers have explored the HT effects on the structural responses of the ME composites. Like the conventional composites, the physical and material properties of ME composites also vary when exposed to the HT environment. The additional interaction of the HT fields with the rest of the fields results in pyro-coupling [25–28]. Zhao et al. [29] studied the coupled response of porous functionally graded PEM structures using finite element methods. Akbarzadeh and Chen [30] evaluated the influence of thermal loads on the static response of rotating PEM cylinders. The vibration response of thermally loaded PEM sandwich nanoplates was probed by Koc et al. [31]. Zheng et al. [32] made use of the analytical formulation and assessed the influence of thermal loads on the nonlinear deflections of ME plates. Zhou and Qu [33] incorporated NURBS technique and studied the static and dynamic response of ME beams. Zhou et al. [34] and Nie et al. [35] investigated the effect of HT loads on the variation of static parameters and dynamic parameters on ME structures using an interpolation method. Ni et al. [36] probed the buckling characteristics of hygrothermally loaded ME cylindrical structures. Using modified coupled stress theory, Zheng et al. [37] examined the nonlinear post-buckling behavior of ME microbeams exposed to a thermal environment. Based on a similar approach, Gui and Wu [38] investigated the buckling response of ME nano cylindrical shells operated in the thermal environment. Esen and Ozmen [39] investigated the influence of thermal fields on the frequency response of ME nanoplates using nonlocal strain gradient theory (NSGT). Sirimontree et al. [40] investigated the sound transmission behavior of sandwich ME composites using NSGT. Barati and Shariyat [41] developed an analytical model to determine the static performance of ME composite structures exposed to HT loads. The coupled vibration response of ME composites with CNTs exposed to HT environment was analytically assessed by Dat et al. [42].

In many instances, various engineering applications demand the structure to be tilted unlike the conventional rectangular/square configuration. The structures with the tilted edges are known as skewed structures. The structural response of skewed structures is unique as opposed to conventional rectangular structures [43–47]. Also, the numerical approach adopted to investigate these structures vary. Tao and Dai [48] assessed the nonlinear frequency response of skewed porous nanoplates. Kallannavar and Kattimani [49] probed on the attenuated frequencies of skewed 3D printed CNT plates. Noroozi and Malekzadeh [50] investigated the effect of moving loads on the deflections of graphene composite skewed plates using meshfree radial interpolation method. Duan et al. [51] examined the stability response of skewed elastic plates using an analytical approach. Based on the isogeometric approach, Sengar et al. [52] probed on the effect of compressive loading condition on the postbuckled vibration behavior of a skewed sandwich plate with a metal core. Abdollahi et al. [53] investigated the influence of geometrical skewness on the aero-flutter behavior of porous plates. The material properties are considered to be functionally graded. The influence of skewness on the frequency of the graphene composite plate was studied by Kiani and Zur [54]. Naveenkumar et al. [55] examined the synergistic effect of porosity and skewness on the nonlinear transient behavior of functionally graded plates. Vaghefi [56] studied the thermos-elasto-plastic behavior of skew plates based on the meshless approach. Farsadi et al. [57] investigated the nonlinear aero-elastic response of skewed plates made of metal/ceramic combination. The influence of geometrical skewness associated with the topological texture on the coupled frequency of PEM-auxetic plates was studied by Mahesh [58].

Though, the computational approaches discussed above are effective analysis tools for the multiphysics structural problems, the computational efforts involved are extremely high. To this end, the predictive tools prove to be extremely useful [59–64]. Among them, reduced order models based on artificial neural networks (ANN) are gaining attention for its ease and accuracy of prediction. ANN becomes very handy when the problem under consideration has too many parameters affecting the output. Turan et al. [65] adopted am ANN-based prediction strategy and assessed the free vibration and buckling behavior of porous beams. Fallah and Agdham [66] employed a physics informed neural network (PINN) approach to predict the bending and vibration characteristics of porous beams. Truong et al. [67] adopted deep forward neural network approach to optimize the static and vibration response of functionally graded beams. Wang and Zhuang [68] proposed a convolutional neural network to predict nonlocal behavior of flexoelectric structures. The prediction of damage in CNT reinforced composite plate was proposed by Le et al. [69] using deep forward neural network. Based on group method of data handling algorithm, the natural frequency of a cracked beam was predicted by Alijani et al. [70]. Chi et al. [71] developed an algorithm based on deep neural network (DNN) and trained it to predict the stability characteristics of the composite structures. Petrolo and Carrera [72] exploited the benefits of neural networks and studied the accuracy of different shell models. Lihua [73] made use of DNN to predict the buckling behavior of hybrid composite panels subjected to HT environment. Very recently, the PEM structures’ coupled response prediction using ANN was proposed by Mahesh and coresearchers [74, 75].

From the exhaustive literature survey carried out, it was realized that, very limited works have been reported on assessing the coupled free vibration characteristics (i.e., coupled frequency) of MME plates under the presence of geometric skewness and hygrothermal environment. In this direction, this work makes a first attempt to integrate these parameters and investigate the frequency response of MME composite plate. To consider all possibilities of stiffness coefficient variation, four forms of material property relation are assumed. In addition, both negative and positive empirical constants dictating the effect of temperature and moisture gradients on the elastic stiffness are considered. For the first time in literature, an artificial neural network (ANN) based computational approach is adopted to predict the natural frequency of MME plates when exposed to different hygrothermal environments, geometrical skewness, volume fractions and material property relations. Since, the considered structural problem is composed of different material, working and geometrical parameters, the proposed ANN based predicting strategy aids the design engineer in accelerating the simulations without compromising on accuracy.

2. Materials and methods

A multifunctional ME composite which exhibits interaction between electric, magnetic and elastic fields is considered in this study. In specific, multiphase ME composites with different percentage of piezoelectric and piezomagnetic phases are considered (Table 1). The constitutive relation between the various fields of ME composites can be represented as follows [26]: (1) $$\begin{array}{c}\left\{\sigma \right\}\text{=}\left[C\right]\text{{}\epsilon \text{}}-\left[e\right]\left\{E\right\}-\left[q\right]\left\{H\right\}\\ \left\{D\right\}={\left[e\right]}^T\left\{\epsilon \right\}+\left[\eta \right]\left\{E\right\}+\left[m\right]\left\{H\right\}\\ \left\{B\right\}={\left[q\right]}^T\left\{\epsilon \right\}+\left[m\right]\left\{E\right\}+\left[\mu \right]\left\{H\right\}\end{array}$$(1)

Table 1. Material properties corresponding to different volume fraction V_f of BaTiO₃–CoFe₂O₄ [26].

Material property	Material constants	0 Vf	0.2 Vf	0.4 Vf	0.5 Vf	0.6 Vf	0.8 Vf	1 Vf
Elastic constants (GPa)	C11=C22	286	250	225	220	200	175	166
	C12	173	146	125	120	110	100	77
	C13=C23	170	145	125	120	110	100	78
	C33	269.5	240	220	215	190	170	162
	C44=C55	45.3	45	45	45	45	50	43
	C66	56.5	52	50	50	45	37.5	44.5
Piezoelectric constants (C/m^2)	e31	0	−2	−3	−3.5	−3.5	−4	−4.4
	e33	0	4	7	9.0	11	14	18.6
	e15	0	0	0	0	0	0	11.6
Dielectric constant (10^−9 C^2/Nm^2)	ε11=ε22	0.08	0.33	0.8	0.85	0.9	1	11.2
	ε33	0.093	2.5	5	6.3	7.5	10	12.6
Magnetic permeability (10^−4 Ns^2/C^2)	μ11=μ22	−5.9	−3.9	−2.5	−2.0	−1.5	−0.8	0.05
	μ33	1.57	1.33	1	0.9	0.75	0.5	0.1
Piezomagnetic constants (N/Am)	q31	580	410	300	350	200	100	0
	q33	700	550	380	320	260	120	0
	q15	560	340	220	200	180	80	0
Magneto-electric constant (10^−12Ns/VC)	m11=m22	0	2.8	4.8	5.5	6	6.8	0
	m33	0	2000	2750	2600	2500	1500	0
Density (kg/m^3)	ρ	5300	5400	5500	5550	5600	5700	5800

Meanwhile, in the current work, the elastic stiffness coefficients are assumed to be dependent on the thermal and moisture fields. The influence of temperature and moisture gradients are related to the elastic stiffness coefficients using empirical constants γ and β, respectively. Henceforth, the following four types of elastic stiffness material property (MP) relations are used to define the variation of elastic stiffness coefficients, (2a) $$\text{MP-1:}\left[C\right]=\left[C_0\right]\left(1+\gamma \Delta T+\beta \Delta m\right)$$(2a) (2b) $$\text{MP-2:}\left[C\right]=\left[C_0\right]\left(1-\gamma \Delta T-\beta \Delta m\right)$$(2b) (2c) $$\text{MP-3:}\left[C\right]=\left[C_0\right]\left(1+\gamma \Delta T-\beta \Delta m\right)$$(2c) (2d) $$\text{MP-4:}\left[C\right]=\left[C_0\right]\left(1-\gamma \Delta T+\beta \Delta m\right)$$(2d)

The schematics of the ME plate is shown in Figure 1a and its plate kinematics is based on Reddy’s higher order shear deformation theory [76]. The components of displacement at any point on the plate is related to the mid-plane displacements and rotations as follows: (3) $$\begin{array}{c}u=u_0+z{\theta }_x-\frac{4}{3h^2}{z}^3\left({\theta }_x+\frac{\partial w_0}{\partial x}\right)\\ v=v_0+z{\theta }_y-\frac{4}{3h^2}{z}^3\left({\theta }_y+\frac{\partial w_0}{\partial y}\right)\\ w=w_0\end{array}$$(3)

Figure 1. Schematic of (a) ME composite plate subjected to hygrothermal environment (b) geometrical skewness.

The governing equations of motion to determine the natural frequencies of ME composite is derived using Hamilton’s principle. For a ME composite, the Hamilton’s principle can be expressed using kinetic and potential energies as [28], (4) $$\begin{array}{c}\delta T_p=\frac{1}{2}\left({\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_b\right\}}^T\left\{{\sigma }_b\right\}+{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_s\right\}}^T\left\{{\sigma }_s\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{E\right\}}^T\left\{D\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{H\right\}}^T\left\{B\right\}\right)d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_k={\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta \left\{d_t\right\}\rho \left\{\ddot{d_t}\right\}d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_p+\delta T_k=0\end{array}$$(4)

Further, the constitutive relations, strains, electric displacement, and magnetic flux density appearing in Eq. (4)(4) $$\begin{array}{c}\delta T_p=\frac{1}{2}\left({\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_b\right\}}^T\left\{{\sigma }_b\right\}+{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_s\right\}}^T\left\{{\sigma }_s\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{E\right\}}^T\left\{D\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{H\right\}}^T\left\{B\right\}\right)d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_k={\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta \left\{d_t\right\}\rho \left\{\ddot{d_t}\right\}d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_p+\delta T_k=0\end{array}$$(4) are expressed using the FE entities. In this regard, the plate is discretized using eight noded isoparametric elements with nine degrees of freedom, which can be listed as three translational ($\left\{d_{ti}\right\}={\left[u_i{\text{ v}}_i{w}_i\right]}^T$), two rotational ($\left\{d_{ri}\right\}={\left[{\theta }_x{\theta }_y\right]}^T$), two higher order rotational ($\left\{d_{r*i}\right\}={\left[{\nu }_x{\nu }_y\right]}^T$) and one for each electric potential ($\varphi$) and magnetic potential ($\psi$) [77]. Meanwhile, the bending and shear components of strains can be represented as [77]: (5) $$\begin{array}{c}\left\{{\epsilon }_b\right\}=\left[B_{tb}\right]\left\{d_t^e\right\}+z\left[B_{rb}\right]\left\{d_r^e\right\}+c_1{z}^3\left[B_{rb}\right]\left\{d_r^e\right\}+c_1{z}^3\left[B_{rb}\right]\left\{d_{r\ast }^e\right\}\\ \left\{{\epsilon }_s\right\}=\left[B_{ts}\right]\left\{d_t^e\right\}+\left[B_{rs}\right]\left\{d_r^e\right\}+c_2{z}^2\left[B_{rs}\right]\left\{d_r^e\right\}+c_2{z}^2\left[B_{rs}\right]\left\{d_{r\ast }^e\right\}\end{array}$$(5)

The different strain-displacement matrices of Eq. (5) are represented in the Appendix. Considering quasi-static behavior of magnetic and electric waves, the relation between electric intensity and magnetic intensity (H) to potential gradients [$\varphi$ and $\psi ],$ can be expressed as follows [44]: (6) $$\begin{array}{c}\left\{E\right\}={\left\{\begin{array}{ccc}{E}_x& E_y& E_z\end{array}\right\}}^T=\left[B_{\varphi }\right]\left\{{\varphi }^e\right\}\\ \left\{H\right\}={\left\{\begin{array}{ccc}{H}_x& H_y& H_z\end{array}\right\}}^T=\left[B_{\psi }\right]\left\{{\psi }^e\right\}\end{array}$$(6)

Substituting Eqs. (1)(1) $$\begin{array}{c}\left\{\sigma \right\}\text{=}\left[C\right]\text{{}\epsilon \text{}}-\left[e\right]\left\{E\right\}-\left[q\right]\left\{H\right\}\\ \left\{D\right\}={\left[e\right]}^T\left\{\epsilon \right\}+\left[\eta \right]\left\{E\right\}+\left[m\right]\left\{H\right\}\\ \left\{B\right\}={\left[q\right]}^T\left\{\epsilon \right\}+\left[m\right]\left\{E\right\}+\left[\mu \right]\left\{H\right\}\end{array}$$(1) , (2), (5)(5) $$\begin{array}{c}\left\{{\epsilon }_b\right\}=\left[B_{tb}\right]\left\{d_t^e\right\}+z\left[B_{rb}\right]\left\{d_r^e\right\}+c_1{z}^3\left[B_{rb}\right]\left\{d_r^e\right\}+c_1{z}^3\left[B_{rb}\right]\left\{d_{r\ast }^e\right\}\\ \left\{{\epsilon }_s\right\}=\left[B_{ts}\right]\left\{d_t^e\right\}+\left[B_{rs}\right]\left\{d_r^e\right\}+c_2{z}^2\left[B_{rs}\right]\left\{d_r^e\right\}+c_2{z}^2\left[B_{rs}\right]\left\{d_{r\ast }^e\right\}\end{array}$$(5) and (6)(6) $$\begin{array}{c}\left\{E\right\}={\left\{\begin{array}{ccc}{E}_x& E_y& E_z\end{array}\right\}}^T=\left[B_{\varphi }\right]\left\{{\varphi }^e\right\}\\ \left\{H\right\}={\left\{\begin{array}{ccc}{H}_x& H_y& H_z\end{array}\right\}}^T=\left[B_{\psi }\right]\left\{{\psi }^e\right\}\end{array}$$(6) in Eq. (4)(4) $$\begin{array}{c}\delta T_p=\frac{1}{2}\left({\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_b\right\}}^T\left\{{\sigma }_b\right\}+{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{{\epsilon }_s\right\}}^T\left\{{\sigma }_s\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{E\right\}}^T\left\{D\right\}-{\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta {\left\{H\right\}}^T\left\{B\right\}\right)d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_k={\int }_{{\mathrm{\Omega }}^{\mathit{plate}}}^{}\delta \left\{d_t\right\}\rho \left\{\ddot{d_t}\right\}d{\mathrm{\Omega }}^{\mathit{plate}}\\ \delta T_p+\delta T_k=0\end{array}$$(4) , and performing condensation approach the equations of motion dictating the free vibrations of ME composite plate can be written as [44]: (7) $$\begin{array}{c}\left[M_{tt}^e\right]\left\{{\ddot{d}}_t^e\right\}+\left[K_{tt}^e\right]\left\{d_t^e\right\}+\left[K_{tr}^e\right]\left\{d_r^e\right\}+\left[K_{t{r}^{\mathrm{*}}}^e\right]\left\{d_{r^{\mathrm{*}}}^e\right\}+\left[K_{t\varphi }^e\right]\left\{{\varphi }^e\right\}+\left[K_{t\psi }^e\right]\left\{{\psi }^e\right\}=0\\ \left[K_{rt}^e\right]\left\{d_t^e\right\}+\left[K_{rr}^e\right]\left\{d_r^e\right\}+\left[K_{tr}^e\right]\left\{d_{r^{\mathrm{*}}}^e\right\}+\left[K_{r\varphi }^e\right]\left\{{\varphi }^e\right\}+\left[K_{r\psi }^e\right]\left\{{\psi }^e\right\}=0\\ \begin{array}{c}\left[K_{{rt}^{\mathrm{*}}}^e\right]\left\{d_t^e\right\}+\left[K_{r{r}^{\mathrm{*}}}^e\right]\left\{d_r^e\right\}+\left[\left[K_{r^{\mathrm{*}}{r}^{\mathrm{*}}}^e\right]\right]\left\{d_{r^{\mathrm{*}}}^e\right\}+\left[K_{r^{\mathrm{*}}\varphi }^e\right]\left\{{\varphi }^e\right\}+\left[K_{r^{\mathrm{*}}\psi }^e\right]\left\{{\psi }^e\right\}=0\\ \left[K_{\varphi t}^e\right]\left\{d_t^e\right\}+\left[K_{r\varphi }^e\right]\left\{d_r^e\right\}+\left[K_{{\varphi r}^{\mathrm{*}}}^e\right]\left\{d_{r^{\mathrm{*}}}^e\right\}-\left[K_{\varphi \varphi }^e\right]\left\{{\varphi }^e\right\}\left[K_{\varphi \psi }^e\right]\left\{{\psi }^e\right\}=0\\ \left[K_{\psi t}^e\right]\left\{d_t^e\right\}+\left[K_{r\psi }^e\right]\left\{d_r^e\right\}+\left[K_{{\psi r}^{\mathrm{*}}}^e\right]\left\{d_{r^{\mathrm{*}}}^e\right\}-\left[K_{\psi \varphi }^e\right]\left\{{\varphi }^e\right\}\left[K_{\psi \psi }^e\right]\left\{{\psi }^e\right\}=0\end{array}\end{array}$$(7)

Meanwhile, to incorporate the effect of geometrical skewness as shown in Figure 1b, the axes and the degrees of freedom (DoF) along the skewed edges are transformed with respect to the skew angle. Hence, the DoF used in the formulation can now be rewritten using the skew degrees of freedom ($\left\{d_t^1\right\},$ $\left\{d_r^1\right\}$ and $\left\{d_{r*}^1\right\}$) as follows: (8) $$\left\{d_t\right\}=\left[T_t\right]\left\{d_t^1\right\};\left\{d_r\right\}=\left[T_r\right]\left\{d_r^1\right\};\left\{d_{r*}\right\}=\left[T_{r*}\right]\left\{d_{r*}^1\right\}$$(8) here, (9a) $$\left[T_t\right]=\left[\begin{array}{ccc}\chi & \zeta & 0\\ -\zeta & \chi & 0\\ 0& 0& 1\end{array}\right];\left[T_r\right]=\left[\begin{array}{cc}\chi & \zeta \\ -\zeta & \chi \end{array}\right];$$(9a) where (9b) $$\begin{array}{c}\chi =\mathrm{Cos}\phi \\ \zeta =\mathrm{Sin}\phi \end{array}$$(9b)

This transformation in turn alters the corresponding stiffness values and hence the natural frequency. The transformed stiffness matrices relating to different combinations of the DoFs can be rewritten as: (10) $$\begin{array}{c}\left[\overline{K_{tt}^e}\right]={\left[T_1\right]}^T\left[K_{tt}^e\right]\left[T_1\right];\left[\overline{K_{tr}^e}\right]={\left[T_1\right]}^T\left[K_{tr}^e\right]\left[T_2\right];\left[\overline{K_{tr*}^e}\right]={\left[T_1\right]}^T\left[K_{tr*}^e\right]\left[T_2\right];\\ \left[\overline{K_{rr}^e}\right]={\left[T_2\right]}^T\left[K_{rr}^e\right]\left[T_2\right];\left[\overline{K_{rr*}^e}\right]={\left[T_2\right]}^T\left[K_{rr*}^e\right]\left[T_2\right];\left[\overline{K_{r*r*}^e}\right]={\left[T_2\right]}^T\left[K_{r*r*}^e\right]\left[T_2\right];\\ \left[\overline{M_{}^e}\right]={\left[T_1\right]}^T\left[M_{}^e\right]\left[T_1\right];\end{array}$$(10) where $\left[T_1\right]$ and $\left[T_2\right]$ are given by (11) $$\begin{array}{c}\left[T_1\right]=\left[\begin{array}{cccccccc}\left[T_t\right]& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}\\ \overline{o}& \left[T_t\right]& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}\\ \overline{o}& \overline{o}& \left[T_t\right]& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}\\ \overline{o}& \overline{o}& \overline{o}& \left[T_t\right]& \overline{o}& \overline{o}& \overline{o}& \overline{o}\\ \overline{o}& \overline{o}& \overline{o}& \overline{o}& \left[T_t\right]& \overline{o}& \overline{o}& \overline{o}\\ \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \left[T_t\right]& \overline{o}& \overline{o}\\ \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \left[T_t\right]& \overline{o}\\ \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \overline{o}& \left[T_t\right]\end{array}\right]\\ \left[T_2\right]=\left[\begin{array}{cccccccc}\left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]& \overline{\overline{o}}\\ \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \overline{\overline{o}}& \left[T_r\right]\end{array}\right]\end{array}$$(11) where $\overline{o}$ and $\overline{\overline{o}}$ are (3 × 3) and (2 × 2) identity matrices, respectively.

Once the transformations have been carried out, the coupled equilibrium equation of motion for free vibration analysis can be explicitly represented as follows: (12) $$\left[M_{tt}\right]\left\{{\ddot{d}}_t\right\}+\left[K_{eq}\right]\left\{d_t\right\}=0$$(12) where $\left[K_{eq}\right]$ is the equivalent stiffness matrix.

The detailed condensation procedure and explicit expressions of the various stiffness and rigidity matrices contributing to Eq. (12)(12) $$\left[M_{tt}\right]\left\{{\ddot{d}}_t\right\}+\left[K_{eq}\right]\left\{d_t\right\}=0$$(12) is incorporated from the authors’ own work [44].

3. Artificial neural network (ANN) model

ANN mimics the human neuron system which makes decision or prediction with the aid of weights, bias and an activation function, whose basic structure is shown in Figure 2a. In the preliminary step, the weight functions (w_ij) are multiplied with the corresponding inputs (y_i) and its summation is added to a non-zero bias (b_i) to obtain S_i, as follows: (13) $$S_i=\sum _{i=1}^n{w}_{ij}{y}_i+b_i$$(13)

Figure 2. The schematics of (a) basic elements constituting an ANN model (b, c) architecture of the ANN model developed to predict vibration characteristics of SMEE plates.

Further, the output parameter z_i can be obtained by enforcing a transfer or activation function f on S_i, as follows: (14) $$z_i=f\left(S_i\right)$$(14)

Meanwhile, the statistical performance of the ANN can be assessed through mean square error (MSE) and the coefficient of determination which can be mathematically represented as follows: (15) $$\begin{array}{c}MSE=\left[\frac{1}{n}\sum _{i=1}^n{\left(Y_i-Y_{id}\right)}^2\right]\\ R=\sqrt{1-\frac{\sum _{i=1}^n{\left(Y_i-Y_{id}\right)}^2}{\sum _{i=1}^n{\left(Y_{id}-Y_m\right)}^2}}\end{array}$$(15) here, Y_id, Y_i, Y_m denote the actual, predicted, and average values, respectively. Further, n denotes the number of points.

In the present research, the ANN model has been created in MATLAB software. To train the ANN model, more than a million data points related to eight inputs (Table 2) collected from the FE simulations were used. The total time taken to simulate the data using FEM was 1258 min. The dataset is split into 70, 15 and 15% for training, validation and testing, respectively. The Levenberg–Marquardt algorithm is used in this current study. Meanwhile, the most optimum number of neurons was decided based on the trial-and-error method. The schematic representation of the ANN model with input, hidden and output layers is shown in Figure 2b.

Table 2. Different parameters and their ranges selected for the development of ANN model.

Sl. No	Parameters	Range
1	Volume fraction (Vf)	0, 20, 40, 60, 80, 100%
2	Skew angle	0, 15, 30, 45, 60
3	Empirical constant (γ, β)	± 0.0025, ± 0.005, ± 0.01
4	Material property relation	MP-1, MP-2, MP-3, MP-4
5	Temperature and moisture gradients	ΔT = 0–100 (in intervals of 10) Δm = 0–2% (in intervals of 0.05%)
6	Boundary condition	CCCC, SSSS, FSFS, CFFF, CCFF, FCFC, SFSF, SSSF, SCSF, SCSS, SCSC, CFCF, CCCF
7	a/h ratio	10-100 (in intervals of 10)
8	a/b ratio	0.5-3 (in intervals of 0.5)

4. Results and discussion

The results of the ANN model together with FE results are presented in this section. The regression plots corresponding to different neuron numbers are shown in Figure 3a. From this figure and Table 3, the optimum number of neurons was selected to be thirty considering the R and RMSE values. Figure 3b depicts the regression results obtained from the trained ANN model with thirty neurons. It can be noticed that a good correlation is achieved. Figure 4a, b depict the error histogram and the performance of the ANN model developed.

Figure 3. (a) Comparison of the neuron number on the overall value of ‘R’ of the ANN model (b) regression plots for training, and testing of the optimum ANN model with thirty neurons and Levenberg-Marquardt algorithm.

Figure 4. The plot of (a) Error histogram (b) performance of the ANN model developed to predict vibration characteristics of ME plates.

Table 3. The performance parameters of the ANN model.

Parameters
No of neurons	5	10	15	20	25	30
Computational time (in minutes)	3.28	11.32	18.45	31.22	72	138
RMSE (training)	10.22%	7.81%	4.75%	3.23%	2.77%	1.95%
RMSE (testing)	15.43%	12.17%	7.36%	4.89%	3.74%	2.3%
R	0.9929	0.9982	0.9992	0.9994	0.9994	0.9998
Epoch	423	611	935	1000	1000	1000

The variation of the coupled fundamental frequency of MME composite skew plate operated in different hygrothermal environments and possessing various elastic stiffness material property (MP) relations has been investigated in this section. In addition, a comparison between the FE simulated and ANN predicted results has been depicted to demonstrate the credibility of both the approaches proposed in this work. Firstly, the convergence study followed by validation is presented. Also, from Tables 3 and 4 it is affirmed that the natural frequencies of ME rectangular plate with alternate piezoelectric and piezomagnetic layers as reported in Refs. [78, 79] closely agree with the FE simulated results with converged 10 × 10 mesh size. Henceforth, the proposed FE model can be used to assess the influence of hygrothermal environment on MME composite plates. Unless and until specified, the following parameters are assumed for simulation: a/h = 100; a/b = 1; clamped (CCCC) boundary condition; skew angle φ = 0^ο; ΔT = 25 K; Δm = 2%; γ = β = 0.005; V_f = 60%

Table 4. Convergence of natural frequency of MEE plate with mesh size.

Mode No.	Mesh size				Moita et al. [78]
	4 × 4	6 × 6	8 × 8	10 × 10
1	2422.14	2428.60	2432.14	2432.21	2449.7
2	6255.74	6264.67	6268.31	6268.40	6280.6
3	6255.74	6264.67	6268.31	6268.40	6280.6

The MME composites is constituted of different volume fractions of piezoelectric and piezo magnetic materials which directly affects the coupled material properties of MME composites. Further, from Eqs. 2(a)–(d)(2a) $$\text{MP-1:}\left[C\right]=\left[C_0\right]\left(1+\gamma \Delta T+\beta \Delta m\right)$$(2a) , it is evident that the elastic stiffness coefficients of the MME composites are directly dependent on the empirical constants and form of elastic stiffness MP. This, in turn, have a significant influence on the natural frequency of the MME composites plate. To this end, an attempt has been made to assess the variation of fundamental natural frequency with changing volume fraction, empirical constants, and MP equation. From Figure 5, it can be seen that with higher volume fraction of ME composite, the natural frequency significantly reduces, which is due to the higher volume fraction of piezoelectric phase with a lower elastic stiffness coefficient. In addition, it is revealed that MP-1 and MP-2 has a predominant effect on the fundamental natural frequency when the empirical constants are positive and negative, respectively. This is because, in the corresponding cases, the corresponding empirical constants result in enhanced elastic stiffness coefficients in accordance with Eqs. 2a(2a) $$\text{MP-1:}\left[C\right]=\left[C_0\right]\left(1+\gamma \Delta T+\beta \Delta m\right)$$(2a) and b. For all the cases, the ANN predicted values closely agree with the FEM simulated results.

Figure 5. Influence of volume fraction and elastic stiffness material property relation on the fundamental frequency of ME plate subjected to hygrothermal environment.

The effect of geometrical skewness on the natural frequency of the ME plate is analyzed in Figure 6a considering temperature and moisture independent MP. In order to evaluate the generality of the trained ANN model, the skew angles not used for training are considered. Due to the reduction in mass of the plate with higher skew angle, the natural frequency tends to increase. Moreover, it is evident that the ANN predicted values are in close agreement with the FEM simulated results. Also, the effect of volume fraction appears to be predominant when the plate has a higher skew angle. This may be attributed to the coupling effects on the plate’s stiffness. Further, considering the ME plate with temperature and moisture dependent elastic stiffness coefficients, the influence of skew angle is probed. For the considered values of γ and β, the variation of fundamental frequency with geometrical skewness is plotted in Figure 6b. It can be seen that the geometrical skewness has a significant effect when the plate’s elastic stiffnesses follow MP-1 type of relation. This is because for MP-1 relation, both the thermal and hygral fields contribute to improving the elastic stiffness coefficients.

Figure 6. Influence of (a) volume fraction and (b) skew angle on the fundamental frequency of ME plate.

The influence of positive and negative values of empirical constant ‘γ’ associated with the volume fraction and elastic stiffness MP relation on the fundamental frequency of ME plate is shown in Figures 7 and 8, respectively. It can be inferred from these figures that for positive values of γ, the natural frequencies of ME plate with MP-1 and MP-3 relations increase with the higher values of γ. The trend is reversed for MP-2 and MP-4 relations. Contrarily, for negative values of empirical constant γ, the plate’s frequency improves with lesser negative γ for MP-1 and MP-3 relations but reduces for MP-2 and MP-4 relations. This may be directly related to the variation in the elastic stiffness coefficients. On the other hand, the study is extended to investigate the effect of β on the coupled frequencies. From Figure 9, it can be seen that as opposed to γ, the influence of β is very minimal. This may be because of Eqs. 2(a)–(d)(2a) $$\text{MP-1:}\left[C\right]=\left[C_0\right]\left(1+\gamma \Delta T+\beta \Delta m\right)$$(2a) which dictates the influence of the thermal field over the hygral field. From Figures 7–9, the accuracy of the ANN model to predict the fundamental frequency for different cases of MP equations and empirical constants is also justified (Table 5).

Figure 7. Influence of empirical constant γ (positive) on the fundamental frequency of ME plate with (a) MP-1 (b) MP-2 (c) MP-3 (d) MP-4 material property relations and volume fractions.

Figure 8. Influence of empirical constant γ (negative) on the fundamental frequency of ME plate with (a) MP-1 (b) MP-2 (c) MP-3 (d) MP-4 material property relations and volume fractions.

Figure 9. Influence of empirical constant β (positive) on the fundamental frequency of ME plate with (a) MP-1 (b) MP-2 (c) MP-3 (d) MP-4 material property relations and volume fractions.

Table 5. Verification of the present FE formulation for natural frequencies of layered ME plates (h = 0.3 m; a = b = 1m).

Stacking sequence	Mode No.	Non-dimensional frequency		% Error
		Chen et al. [79]	Present
BFB	1	1.3434	1.3482	0.3601
	2	2.2199	2.2279	0.3597
	3	2.2199	2.2279	0.3597
FBF	1	1.4463	1.4346	−0.8109
	2	2.3602	2.3390	−0.8979
	3	2.3602	2.3392	−0.8905

The investigation is extended to evaluate the effect of different magnitudes of thermal and moisture gradient associated with the empirical constants and the material property relations on the fundamental frequency of ME plate. From Figures 10 and 11, it is noticed that for ME plate with positive empirical constants, the higher value of thermal and moisture gradient enriches the frequency for MP-1 and MP-3 type of elastic stiffness coefficient whereas, the frequency linearly decreases for MP-2 and MP-4 type of relation. This entire trend is reversed for negative empirical constants due to the obvious reasons mentioned previously. Meanwhile, it can also be seen that the thermal gradient exhibits a superior effect than the moisture gradient in both cases.

Figure 10. Influence of (a) thermal (b) moisture gradient on the fundamental frequency of ME plate with different material property relation and positive empirical constants.

Figure 11. Influence of (a) thermal (b) moisture gradient on the fundamental frequency of ME plate with different material property relation and negative empirical constants.

From the above analysis, it is evident that the volume fraction, material property relation, geometric skewness and empirical constants significantly affect the coupled stiffness and hence the natural frequency of the ME plate. Hence, it is extremely important to assess the integrated effects of these parameters on the fundamental frequency. Simultaneously considering the influence of volume fraction, material property relation and skew angle, the variation of fundamental frequency against the thermal gradient is plotted in Figure 12. It is easily visualized that the higher volume fraction and lower skew angle have detrimental effects on the fundamental frequency of ME plate. Similarly, from Figure 13, it can be witnessed that the empirical constants have a predominant effect associated with MP-1 relation and lower volume fraction.

Figure 12. Integrated effects of volume fraction, material property relation and geometric skewness on the fundamental frequency of ME plate with (a) φ = 0^ο (b) φ = 15^ο (c) φ = 30^ο (Δm = 2%; CCCC; γ = β = 0.05).

Figure 13. Integrated effects of volume fraction, different empirical constants and geometric skewness on the fundamental frequency of ME plate with (a) MP-1 (b) MP-2 (c) MP-3 (d) MP-4 material property relation (ΔT = 25; Δm = 2%; CCCC).

In all the parametric studies involving geometric and material parameters, the predictions of the ANN model match very well with that of the FE results with a maximum error of 2.3%, including for the plate configurations not seen by the ANN model during training, thus demonstrating its generality.

5. Conclusions

In the present research, an ANN-based machine learning approach is developed to investigate the frequency response of skewed multiphase ME composite plate subjected to an external hygrothermal environment. More than a million data points obtained from in-house developed FE models are used to train an ANN model using Levenberg–Marquardt algorithm. Four different elastic stiffness material property relations have been considered. The temperature and moisture gradients are related to the elastic stiffness using the empirical constants. The proposed FE model used to generate the dataset displays its versatility over the conventional analytical models to simulate frequency of skewed plates. The numerical simulations suggest that the higher volume fraction of MME composites and lower skew angle are detrimental to the frequency of MME composite plate. Alongside, the material property relation and the sign of the empirical constant play a prominent role is deciding the coupled frequency response. Despite its versatility, the FEM simulations are time-consuming owing to various field interactions resulting from the material, geometrical and hygrothermal gradient parameters. In this context, the trained ANN model is found to simulate the results 6.34 times quicker than the FEM tool. It is believed that with more complex operating conditions and geometrical configurations, the computational efficiency and benefits of the ANN model enhances. Hence, the use of ANN-based reduced-order models can be an appropriate approach for structural design and optimization problems involving such multifunctional composite materials with multiple design parameters.

Data availability

The raw/processed data required to reproduce these findings cannot be shared at this time as the data also forms part of an ongoing study.

Additional information

Funding

The financial support by The Royal Society of London through Newton International Fellowship (NIF\R1\212432) is sincerely acknowledged by the authors Vinyas Mahesh and Sathiskumar Anusuya Ponnusami.

(A1)

Appendix

$$\begin{array}{c}{\nu }_x=\frac{\partial w_0}{\partial x};{\nu }_y=\frac{\partial w_0}{\partial y};c_1=-\frac{4}{3h^2};\left\{d_t^e\right\}=\left\{\begin{array}{c}{u}_0\\ v_0\\ w_0\end{array}\right\};\left\{d_r^e\right\}=\left\{\begin{array}{c}{\theta }_x\\ {\theta }_y\end{array}\right\};\left\{d_{r\ast }^e\right\}=\left\{\begin{array}{c}{\nu }_x\\ {\nu }_y\end{array}\right\}\\ \left[B_{tb}\right]=\left[\begin{array}{ccc}\frac{\partial }{\partial x}& 0& 0\\ 0& \frac{\partial }{\partial y}& 0\\ \frac{\partial }{\partial y}& \frac{\partial }{\partial x}& 0\end{array}\right];\left[B_{rb}\right]=\left[\begin{array}{cc}\frac{\partial }{\partial x}& 0\\ 0& \frac{\partial }{\partial y}\\ \frac{\partial }{\partial y}& \frac{\partial }{\partial x}\end{array}\right]\\ c_2=-\frac{4}{h^2};\left[B_{ts}\right]=\left[\begin{array}{ccc}0& 0& \frac{\partial }{\partial x}\\ 0& 0& \frac{\partial }{\partial y}\end{array}\right];\left[B_{rs}\right]=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right]\end{array}$$