Improved Finite Element Thermomechanical Analysis of Laminated Composite and Sandwich Plates Using the New Enhanced First-Order Shear Deformation Theory

: This paper proposes a simple yet accurate finite element (FE) formulation for the ther-momechanical analysis of laminated composites and sandwich plates. To this end, an enhanced first-order shear deformation theory including the transverse normal effect based on the mixed variational theorem (EFSDTM_TN) was employed in the FE implementation. The primary objective of the FE formulation was to systematically interconnect the displacement and transverse stress fields using the mixed variational theorem (MVT). In the MVT, the transverse stress field is derived from the efficient higher-order plate theory including the transverse normal effect (EHOPT_TN), to enhance the solution accuracy, whereas the displacement field is defined by the first-order shear deformation theory including the transverse normal effect (FSDT_TN), to amplify the numerical efficiency. Furthermore, the transverse displacement field is modified by incorporating the components of the external temperature loading, enabling the consideration of the transverse normal strain effect without introducing additional unknown variables. Based on the predefined relationships, the proposed FE formulation can extract the C0-based computational benefits of FSDT_TN, while improving the solution accuracy for thermomechanical analysis. The numerical performance of the proposed FE formulation was demonstrated by comparing the obtained solutions with those available in the literature, including 3-D exact solutions.


Introduction
In recent years, the utilization of high-strength and lightweight structures has continued to improve energy efficiency in line with a wide range of environmental issues.In this regard, fiber-reinforced composite materials capable of providing an optimized stiffness-toweight ratio through a synergistic combination of two or more materials, such as reinforcing fibers and resins, are attracting considerable attention as prospective next-generation materials in various engineering fields.Continuous fiber-reinforced composites are widely employed in various high-value industries, including automotive, civil, and aerospace, owing to their ability to achieve excellent structural properties and multifunctional characteristics.Despite the aforementioned advantages, the distribution of transverse stress in laminated composite structures can give rise to inherent mechanical defects, notably layer slip and delamination.Therefore, accurate prediction of transverse stress is a crucial concern in the structural design process of laminated composite structures [1,2].
Over the past half-century, a range of analysis models based on the equivalent singlelayer theory have been developed to precisely elucidate the transverse behaviors of laminated composite plates .Starting from the well-known classical laminated plate theory (CLPT) and progressing through the first-order shear deformation theory (FSDT), a series of higher-order polynomial theories, including the higher-order shear deformation theory (HSDT), have been developed sequentially [3][4][5][6][7][8].However, most of these theories exhibit limitations in predicting interlaminar stresses because of their inability to enforce transverse shear stress conditions at both the surface and layer interfaces.To address this issue, a series of refined zigzag theories (EHOPT: efficient higher order plate theory, RHSDT: refined higher-order shear deformation theory, RZT: refined zigzag theory) have been proposed [14][15][16][17][18].These theories yield reliable results in predicting the global and local behaviors of laminated composites and sandwich structures by introducing a zigzag displacement field that varies discontinuously at the interlaminar interfaces.However, it requires the use of a nonconventional C1-class shape function (a slope continuity condition along the boundary of the element) in the finite element (FE) formulation, which is incompatible with commercial FE software such as ANSYS 2023 R1 (Ansys, PA, USA) and ABAQUS 2022 (Dassault Systemes, Pairs, France).As an attractive scheme to circumvent C1-class problems in FE analysis, enhanced analysis models (EFSDT: enhanced first-order shear deformation theory, EFSDTM: enhanced first-order shear deformation theory based on mixed variational theorem) were developed for the multiphysics analysis of laminated composites and sandwich plates [19][20][21].Enhanced theories can simultaneously improve the solution accuracy and computational efficiency by systematically deriving reasonable energy relationships between the conventional FSDT and EHOPT.Consequently, these theories allow for a C0-based finite element formulation based on an FSDT-like governing equation, providing clear advantages in terms of compatibility with commercial finite element (FE) software (ANSYS 2023 R1 andABAQUS 2022).
With technological advancements, laminated composites and sandwich structures can be exposed to various external environments, and there is a need to predict their multiphysical behaviors during the design process.In high-temperature environments, thermal deformation and stress induce significant defects.Consequently, thermomechanical analysis should be considered to ensure reliable design solutions for laminated composites and sandwich structures exposed to such conditions.Transverse normal deformation is a very important consideration in thermal analysis.Therefore, well-known analysis models (CLPT, FSDT, HSDT, EHOPT, etc.) that ignore the transverse normal strain effect are not suitable for predicting the thermal behavior of laminated composites and sandwich structures.In this regard, many refined theories have been proposed to precisely describe the thermomechanical responses of laminated composites and sandwich structures .As a higher-order polynomial model, the Lo-Christensen-Wu (LCW) theory attempts to consider the transverse normal strain effect effectively by assuming a smooth parabolic form of the transverse displacement field [29].Furthermore, various refined higher-order and zigzag theories have been proposed for the thermomechanical analysis of laminated composites and sandwich structures [29][30][31][32][33][34][35][36][37][38][39][40][41][42][43][44][45][46][47][48].As one of the most attractive zigzag theories, the efficient higher-order zigzag theory (EHOZT) proposed by Oh and Cho can provide reliable solutions for fully coupled electro-thermo-mechanical problems by enforcing transverse shear stress conditions at both the surface and layer interfaces [40][41][42][43].Kapuria and Achary developed a computationally efficient zigzag theory to predict the thermal behavior of laminated composite structures [44].Although this theory considers the transverse normal strain effect without introducing additional variables into the displacement fields, its applicability in analyzing sandwich plates is limited.This is because the use of different thermal expansion coefficients in adjacent layers can potentially violate the transverse displacement continuity conditions.Among the various enhanced theories for the thermomechanical analysis of laminated composites and sandwich structures [49][50][51], Han et al. proposed an enhanced first-order shear deformation theory including the transverse normal strain effect based on the mixed variational theorem (EFSDTM_TN) to take computational benefits of conventional FSDT [50].The main contribution of EFSDTM_TN is that it considers the transverse normal strain effect without introducing any additional unknown variables by extending the transverse normal displacement field under the prescribed thermal conditions.Furthermore, the transverse displacement continuity conditions are automatically satisfied in the sandwich plates by introducing layer-wise constants.Consequently, EFS-DTM_TN can provide reliable solutions for analyzing the thermomechanical behaviors of laminated composite and sandwich structures while ensuring the computational benefits of the C0-based 5-DOF element in the FE implementation.
To further extend the applicability of the EFSDTM_TN [50], an FE formulation based on the EFSDTM_TN was proposed and numerically tested.An 8-node serendipity element was utilized in the FE formulation to enhance the computational efficiency in deriving the stress distributions.The primary objective of the proposed FE analysis model is to ensure both the solution accuracy and computational efficiency by systematically blending FSDT_TN and EHOPT_TN based on the mixed variational theorem.Furthermore, the thermal responses of laminated composites and sandwich structures can be described more precisely by improving the transverse displacement field.To demonstrate the numerical performance of the proposed FE analysis model, representative thermal-mechanical problems for 2-D laminated composite and sandwich structures were considered as numerical examples.The accuracy and efficiency of the proposed FE analysis model were compared with other numerical results available in the literature, including 3-D exact solutions [57,58] together with the analytical solution of the EFSDTM_TN [50].enhanced theories for the thermomechanical analysis of laminated composites and sandwich structures [49][50][51], Han et al. proposed an enhanced first-order shear deformation theory including the transverse normal strain effect based on the mixed variational theorem (EFSDTM_TN) to take computational benefits of conventional FSDT [50].The main contribution of EFSDTM_TN is that it considers the transverse normal strain effect without introducing any additional unknown variables by extending the transverse normal displacement field under the prescribed thermal conditions.Furthermore, the transverse displacement continuity conditions are automatically satisfied in the sandwich plates by introducing layer-wise constants.Consequently, EFSDTM_TN can provide reliable solutions for analyzing the thermomechanical behaviors of laminated composite and sandwich structures while ensuring the computational benefits of the C0-based 5-DOF element in the FE implementation.

EFSDTM_TN for the
To further extend the applicability of the EFSDTM_TN [50], an FE formulation based on the EFSDTM_TN was proposed and numerically tested.An 8-node serendipity element was utilized in the FE formulation to enhance the computational efficiency in deriving the stress distributions.The primary objective of the proposed FE analysis model is to ensure both the solution accuracy and computational efficiency by systematically blending FSDT_TN and EHOPT_TN based on the mixed variational theorem.Furthermore, the thermal responses of laminated composites and sandwich structures can be described more precisely by improving the transverse displacement field.To demonstrate the numerical performance of the proposed FE analysis model, representative thermal-mechanical problems for 2-D laminated composite and sandwich structures were considered as numerical examples.The accuracy and efficiency of the proposed FE analysis model were compared with other numerical results available in the literature, including 3-D exact solutions [57,58] together with the analytical solution of the EFSDTM_TN [50].

Mixed Variational Theorem
Laminated composites and sandwich plates were considered as numerical models of thermomechanical problems.The geometric shapes and reference coordinates of the laminated plates are shown in Figure 1   In the EFSDTM_TN, displacement and transverse stress fields are assumed independently, with the aim of enhancing both solution accuracy and computational efficiency.These independent fields can then be systematically interconnected based on the mixed variational theorem (MVT).The first variation of the 2-D Hellinger-Reissner functional is defined by Equation (1). where In Equation (1), Ω represents the reference plane of the laminated plates, and mechanical loading ( p i ) is applied to the boundary area (S σ ).Additionally, (•) and (•) * are components defined by the displacement and transverse stress fields, respectively.The mixed part of the MVT plays a critical role in defining reasonable relationships between two independent fields [20,50,51].

Improvement of Transverse Displacement Field
In contrast to the mechanical behavior, the transverse normal strain effect is dominant in the thermal deformation of the laminated composite and sandwich plates.Therefore, this effect should be considered to provide a reliable solution for predicting thermal behavior.Intuitively, assuming a smoothly varying parabolic form for the transverse displacement field can help in this regard.Although this approach is able to predict the thermal behavior of laminated composites and sandwich plates precisely, it involves additional unknown variables as well as complicated 3-D governing equations.Therefore, to provide simple yet accurate solutions for thermal problems, a modified transverse displacement field is introduced as follows [50].
The underlined expressions are newly considered in Equation (3) for a reliable thermal analysis of the laminated composites and sandwich plates.Other than the underlined expressions, this represents a typical transverse displacement field that satisfies the assumption of a plane-stress state (u 3 ≈ u (0) 3 ).Based on Equation (3), the prescribed thermal conditions (T 0 and T 1 ) are utilized to define a smoothly varying parabolic field to consider the transverse normal strain effect without introducing additional unknown variables.Here, N represents the total number of layers and H(x 3 − x 3(k) ) is the Heaviside step function.T 0 and T 1 are the uniform and linear temperature loadings, respectively.In addition, φ (k) is a layer-wise constant that automatically satisfies the plane-stress condition.The value of φ (k) can vary depending on the material composing each layer because it is a function of the material properties and thermal expansion coefficients.Therefore, to fulfill the continuity conditions of u 3 for general layup cases such as sandwich plates, a layer-wise constant (S (k) 3 ) was additionally introduced.This modified form of the transverse displacement field enables simple yet accurate thermomechanical analysis of laminated composites and sandwich plates.

Transverse Stress Field
In this subsection, a reliable transverse stress field based on EHOPT_TN is independently assumed in the MVT to ensure solution accuracy.EHOPT_TN can rigorously satisfy the shear-free conditions at the surface, as well as shear continuity conditions at the layer interfaces by introducing third-order zigzag field in the in-plane displacement field.Furthermore, a modified form of the transverse displacement field was employed to provide reliable solutions to thermomechanical problems.The initial displacement field of EHOPT_TN is expressed as [50] u * α = u * (0) α where α is a linear zigzag field that enforces shear continuity conditions at the layer interfaces.By applying the shear stress conditions to the initial displacement field given in Equation (4), S 3 can be defined through relations in terms of the primary unknown variables (u * (k) i ) and prescribed thermal conditions (T 0 , T 1 ): Detailed definitions of a αβ and b 33 are provided in [50].Furthermore, to satisfy the plate equilibrium state rigorously when applying MVT to the general configuration of laminated structures, in-plane correction factors were introduced in EHOPT_TN.Consequently, Equations ( 4) and ( 5) yield the following refined displacement field for EHOPT_TN [50]: where (φ (N) +φ (1) ) in which where δ αβ is the Kronecker delta function and the in-plane correction factor C N α , C M α can be defined as can be determined by matching the resulting forces and moments in the process of establishing a relationship between the displacement and transverse stress fields based on Saint-Venant's principle [50].Based on the introduction of these in-plane correction factors, it is possible to provide highly reliable solutions for predicting the thermomechanical behavior of laminated composites and sandwich plates.
From Equation ( 6), the transverse stress tensors used in the MVT can be defined as follows:

Displacement Field
Simple displacement field based on FSDT_TN were also considered in MVT to retain computational efficiency [50].The displacement field based on FSDT_TN is given as where the components of FSDT_TN are indicated by overbars to clearly distinguish between the displacement and transverse stress fields in MVT.From Equation ( 11), the strain and in-plane stress tensors used in the MVT can be derived as where α γω and ∆T are the thermal expansion coefficient and temperature distribution.

Relationships between Displacement and Transverse Stress Fields
A reasonable relationship between EHOPT_TN and FSDT_TN can be systematically defined using the mixed part in the MVT as a constraint equation.The related constraint equation is expressed as the following [50]: where γ * 3α and σ * 3α are defined in Equation (10), while γ 3α is defined in Equation ( 12).In the constraint equation, the transverse shear resultant (Q * α ) derived from EHOPT_TN can be expressed as in which Â(0 Equations ( 13)-( 15) yield the relationships between u * (3) α and γ (0) 3α as follows: where Γ (1) in which Consequently, the transverse shear resultant (Q * α ) can be expressed in terms of the FSDT_TN variables by substituting Equation ( 16) into Equation ( 14) as follows: where A * α3β3 , B * α3β3 and D * α3β3 are the effective shear stiffness moduli, which depend on the in-plane correction factors.Thus, these in-plane correction factors and effective shear stiffness moduli should be updated by applying iterative calculations to improve the solution accuracy.Equation (20) indicates that the effective shear correction factor (SCF) can be calibrated automatically using EFSDTM_TN [50].
Based on the reasonable relationship between EHOPT_TN and FSDT_TN, the 2-D Hellinger-Reissner functional can be simplified as Therefore, considering transverse loading ( t 3 ), the governing equations of EFSDTM_TN can be derived as δu and the associated boundary conditions are given by δu It should be noted that the governing equations of EFSDTM_TN are similar to those of conventional FSDT.This implies that the EFSDTM_TN can be extended using a simple FE implementation.
Once the values of all of the unknown variables are determined based on the governing equation, the solution accuracy can be further improved by restoring the displacement field of EHOPT_TN.By applying the least-squares approximation, the following relationships between u * (0) α and u (0) α are obtained [50]: Substituting Equations ( 16) and ( 24) into Equation ( 6), the displacement field of EHOPT_TN can be systematically expressed using only the primary variables of FSDT_TN as follows: αβ Γ (2) in which

Finite Element Formulation Based on EFSDTM_TN
In this section, a finite element formulation based on EFSDTM_TN is presented to further extend its applicability.Considering the stress restoration based on the postprocessing procedure, a well-known 8-node serendipity element was employed in the FE implementation.Based on the FE discretization, the displacement field can be defined by the nodal variables of the 8-node serendipity element as follows [21]: (27) where N i represents the shape function for the (i)th node of the 8-node serendipity element.

Element Stiffness Matrix
The element stiffness matrix can be defined through the principle of minimum potential energy, and all unknown nodal displacements in each 8-node serendipity element can be expressed in vector form as follows [21]: As indicated in Equation (28), each element has 40 degrees of freedom (DOF).
From Equations ( 27)-( 29) and the assumption of small strain-displacement relations, the strain components for each element are defined as where the subscripts (m, b and s) denote the strain components derived from the membrane, bending, and transverse shear parts, respectively.The strain matrices ([B] m , [B] b , and [B] s ) can be written as follows: in which The element stiffness matrix is then reasonably defined by using the aforementioned strain matrices [B] (m,b,s) : in which In Equation ( 34), A , B , D , and G * (0) can be expressed as follows: in which It should be noted that the shear stiffness matrix, G * (0) , is defined based on the effective shear stiffness modulus (A * α3β3 ) instead of G = Q α3β3 .

Extermal Force Vector
In this FE implementation, thermal and mechanical loadings were considered as external force vectors.The corresponding external force vector considering the mechanical loading can be obtained as follows: where and Additionally, based on Equations ( 10) and ( 12), the corresponding external force vector considering thermal loading can be defined as where where T (0,1) , T(0,1) , and α represent vectors consisting of external temperatures and thermal expansion coefficients, respectively, as shown below: It should be also remarked that G * (1) and G * (2) are derived from the effective shear stiffness modulus given in Equation (20).
Based on the above improved FE implementation, both solution accuracy and computational efficiency can be further improved in the process of describing the thermomechanical behaviors of the laminated composite and sandwich plates.

Numerical Results and Discussion
In this section, the numerical performance of the proposed FE analysis model is investigated by considering the characteristic thermomechanical problems of the laminated composites and sandwich plates.For all numerical models, 2-D rectangular laminated plates with simply supported boundary conditions were used as the test beds.The lengthto-thickness ratio of the laminated plates was assumed to be S = L 1 /h = L 2 /h = 4 for mechanical problems and S = L 1 /h = L 2 /h = 5 for thermal problems.
-Each ply of the composite plates for the mechanical problems -Each ply of the composite plates for the thermal problems where the subscripts (•) L and (•) T represent the directions parallel and perpendicular to the fiber configuration.In addition, the material properties of the sandwich plates for thermal and mechanical problems were as follows [50,51]: -Facial sheets of the sandwich plates -Core of the sandwich plates For the thermomechanical problems, the corresponding thermal and mechanical loadings considered when deriving the external force vectors can be expressed as follows: The representative layup configurations of the laminated plates are listed in Table 1.The FE solutions of EFSDTM_TN were then compared with those obtained by conventional C0-class FE analysis models (FSDT, HSDT, and LCW) [3][4][5]8,29] as well as 3-D exact solutions [57,58].The Pagano solutions for thermomechanical problems were considered as benchmark solutions [57,58], and the SCF was assumed to be 5/6 in the conventional FSDT.For reasonable comparison, the numerical results reported herein were normalized in the following form: Numerical results for mechanical problems Numerical results for thermal problems

Validation of the Proposed FE Analysis Model
To examine the numerical errors that may occur during FE analysis, FE solutions based on EFSDTM_TN were validated against those obtained using the analytical approach.To this end, the convergence rate of the FE solutions was numerically verified by comparing the central deflections of the laminated composite and sandwich plates across different mesh densities, as listed in Table 2.The solutions for uniform temperature loading are not compared in Table 2 because of the absence of deflections.
From Table 2, it is observed that the FE solutions gradually converge to the analytical solutions with further refinement of the mesh.In addition, acceptable deflections were obtained when the FE model was discretized into an 8 × 8 mesh density or higher.Although an 8 × 8 mesh density is sufficient to describe the nodal displacement, potential numerical errors could arise when deriving in-plane and transverse shear stresses as these involve higher-order derivatives.The FE solutions for the in-plane and transverse shear stresses of [0 • /90 • /0 • ] laminated composite plates are illustrated in Figure 2. The accuracies of these FE solutions, for various mesh density, were compared with those of analytical solutions.The distributions of the transverse shear stress given in Figure 2 were derived from following 3-D equilibrium equation: analytical solutions.The distributions of the transverse shear stress given in Figure 2 were derived from following 3-D equilibrium equation: Figure 2 shows that an FE model with a 16 × 16 mesh density or higher can yield precise numerical solutions for predicting the local distributions of in-plane and transverse shear stresses.In particular, it can be observed that the FE solutions with a 32 × 32 mesh density can closely approximate the analytical solutions, even for transverse shear stresses that require third-order derivatives.This means that the FE model with a 32 × 32 mesh density can be reasonably applied in thermomechanical analysis of laminated composite and sandwich structures with arbitrary geometry, loading, and boundary conditions.Based on the numerical validation given in Table 2 and Figure 2, the FE solutions for all of the thermomechanical problems were obtained based on a 32 × 32 mesh density to ensure computational accuracy.Considering laminated composites and sandwich structures discretized with a 32 × 32 mesh density, kinematic unknown variables and corresponding total DOFs of the FE models were compared in Table 3. Table 3 shows that the total DOFs of the proposed FE model are the same as FSDT, representing reductions of 45.5% to 55.6% as compared to the total DOFs of the HSDT and LCW, respectively.Therefore, the proposed FE model can clearly improve its computational efficiency in the process of thermo-mechanical analysis.

FE Solutions for the Mechanical Problem
In this subsection, the mechanical behaviors of the laminated composite and sandwich plates are evaluated to verify the numerical performance of the proposed FE model based on EFSDTM_TN.For the mechanical problems, transverse external loading ( 3 0 P ≠ ) was applied to the top surfaces of the laminated plates.Figure 2 shows that an FE model with a 16 × 16 mesh density or higher can yield precise numerical solutions for predicting the local distributions of in-plane and transverse shear stresses.In particular, it can be observed that the FE solutions with a 32 × 32 mesh density can closely approximate the analytical solutions, even for transverse shear stresses that require third-order derivatives.This means that the FE model with a 32 × 32 mesh density can be reasonably applied in thermomechanical analysis of laminated composite and sandwich structures with arbitrary geometry, loading, and boundary conditions.
Based on the numerical validation given in Table 2 and Figure 2, the FE solutions for all of the thermomechanical problems were obtained based on a 32 × 32 mesh density to ensure computational accuracy.Considering laminated composites and sandwich structures discretized with a 32 × 32 mesh density, kinematic unknown variables and corresponding total DOFs of the FE models were compared in Table 3. Table 3 shows that the total DOFs of the proposed FE model are the same as FSDT, representing reductions of 45.5% to 55.6% as compared to the total DOFs of the HSDT and LCW, respectively.Therefore, the proposed FE model can clearly improve its computational efficiency in the process of thermo-mechanical analysis.

FE Solutions for the Mechanical Problem
In this subsection, the mechanical behaviors of the laminated composite and sandwich plates are evaluated to verify the numerical performance of the proposed FE model based on EFSDTM_TN.For the mechanical problems, transverse external loading (P 3 ̸ = 0) was applied to the top surfaces of the laminated plates.
Figure 3 shows the mechanical solutions of the in-plane displacements and stresses for cross-ply laminated composite plates.From Figure 3a, it can be seen that the EFSDTM_TN precisely describes the unsymmetrical zigzag distribution of the in-plane displacement in [0 • /90 • /0 • /90 • ] laminated composite plates.As shown in Figure 3b, EFSDTM_TN can provide a reliable local solution for the in-plane stress by capturing its noncontinuous distribution.However, other theories are only useful for predicting the global behavior of in-plane displacement.4a), it should be noted that EFSDTM_TN can provide a sufficiently reliable solution for the in-plane displacement, even for a sandwich plate.In addition, as shown in Figure 4b, the severe kink distribution of the transverse shear stress is completely captured by EFSDTM_TN.

FE Solutions for the Thermal Problem
To further investigate the numerical capabilities related to the thermal analysis of laminated composites and sandwich plates, several thermal problems were also analyzed.In the case of thermal problems, uniform and linearly distributed temperatures were considered as external loads.For the [0 • /Core/0 • ] sandwich plate, the distributions of the in-plane displacements and transverse shear stresses are shown in Figure 4.In terms of the in-plane displacement (Figure 4a), it should be noted that EFSDTM_TN can provide a sufficiently reliable solution for the in-plane displacement, even for a sandwich plate.In addition, as shown in Figure 4b, the severe kink distribution of the transverse shear stress is completely captured by EFSDTM_TN.4a), it should be noted that EFSDTM_TN can provide a sufficiently reliable solution for the in-plane displacement, even for a sandwich plate.In addition, as shown in Figure 4b, the severe kink distribution of the transverse shear stress is completely captured by EFSDTM_TN.

FE Solutions for the Thermal Problem
To further investigate the numerical capabilities related to the thermal analysis of laminated composites and sandwich plates, several thermal problems were also analyzed.In the case of thermal problems, uniform and linearly distributed temperatures were considered as external loads.

FE Solutions for the Thermal Problem
To further investigate the numerical capabilities related to the thermal analysis of laminated composites and sandwich plates, several thermal problems were also analyzed.In the case of thermal problems, uniform and linearly distributed temperatures were considered as external loads.
The in-plane and transverse shear stresses for a single-layer composite plate under uniform temperature loading (T 0 ̸ = 0) are shown in Figure 5.As stated in Section 2.2, the transverse normal strain effect plays an important role in analyzing the thermal behaviors of laminated composites and sandwich plates.Furthermore, this effect becomes significant under a uniform temperature loading.Considering this aspect, FESDTM_TN and LCW, which reasonably consider the transverse normal strain effect, can accurately describe the thermal stresses of a single-layer composite plate under uniform temperature loading.Another interesting observation from Figure 5 is that the FSDT and HSDT, which cannot consider the transverse normal strain effect, provide meaningless solutions for these thermal stresses.The in-plane and transverse shear stresses for a single-layer composite plate under uniform temperature loading ( 0 0 T ≠ ) are shown in Figure 5.As stated in Section 2.2, the transverse normal strain effect plays an important role in analyzing the thermal behaviors of laminated composites and sandwich plates.Furthermore, this effect becomes significant under a uniform temperature loading.Considering this aspect, FESDTM_TN and LCW, which reasonably consider the transverse normal strain effect, can accurately describe the thermal stresses of a single-layer composite plate under uniform temperature loading.Another interesting observation from Figure 5 is that the FSDT and HSDT, which cannot consider the transverse normal strain effect, provide meaningless solutions for these thermal stresses.Figure 6 compares the thermal distributions of the in-plane displacements and stresses of the cross-ply laminated composite plates under uniform-temperature loading ( 0 0 T ≠ ).As shown in Figure 6, EFSDTM_TN and LCW can precisely capture not only the parabolic distributions of the in-plane displacements, but also the noncontinuous distributions of the in-plane stress, whereas the other theories fail to describe the local distributions of the corresponding thermal behaviors.Figure 6 compares the thermal distributions of the in-plane displacements and stresses of the cross-ply laminated composite plates under uniform-temperature loading (T 0 ̸ = 0).As shown in Figure 6, EFSDTM_TN and LCW can precisely capture not only the parabolic distributions of the in-plane displacements, but also the noncontinuous distributions of the in-plane stress, whereas the other theories fail to describe the local distributions of the corresponding thermal behaviors.The in-plane and transverse shear stresses for a single-layer composite plate under uniform temperature loading ( 0 0 T ≠ ) are shown in Figure 5.As stated in Section 2.2, the transverse normal strain effect plays an important role in analyzing the thermal behaviors of laminated composites and sandwich plates.Furthermore, this effect becomes significant under a uniform temperature loading.Considering this aspect, FESDTM_TN and LCW, which reasonably consider the transverse normal strain effect, can accurately describe the thermal stresses of a single-layer composite plate under uniform temperature loading.Another interesting observation from Figure 5 is that the FSDT and HSDT, which cannot consider the transverse normal strain effect, provide meaningless solutions for these thermal stresses.Figure 6 compares the thermal distributions of the in-plane displacements and stresses of the cross-ply laminated composite plates under uniform-temperature loading ( 0 0 T ≠ ).As shown in Figure 6, EFSDTM_TN and LCW can precisely capture not only the parabolic distributions of the in-plane displacements, but also the noncontinuous distributions of the in-plane stress, whereas the other theories fail to describe the local distributions of the corresponding thermal behaviors.The thermal distributions of the in-plane and transverse shear stresses of the sandwich plates under uniform temperature loading (T 0 ̸ = 0) are shown in Figure 7.As shown in Figure 7a, the thermal solutions obtained by EFSDTM_TN and LCW are in good agreement with the exact solutions by precisely describing the severe noncontinuous distribution of the in-plane stress.In addition, Figure 7b indicates that EFSDTM_TN and LCW provide the best compromised thermal solutions for the local distribution of the transverse shear stress.Consequently, Figures 5-7 demonstrate that the transverse normal strain effect should be considered to accurately describe the thermal behavior of composite materials and sandwich structures subjected to uniform temperature loading.The thermal distributions of the in-plane and transverse shear stresses of the sandwich plates under uniform temperature loading ( 0 0 T ≠ ) are shown in Figure 7.As shown in Figure 7a, the thermal solutions obtained by EFSDTM_TN and LCW are in good agreement with the exact solutions by precisely describing the severe noncontinuous distribution of the in-plane stress.In addition, Figure 7b indicates that EFSDTM_TN and LCW provide the best compromised thermal solutions for the local distribution of the transverse shear stress.Consequently, Figures 5-7 demonstrate that the transverse normal strain effect should be considered to accurately describe the thermal behavior of composite materials and sandwich structures subjected to uniform temperature loading.For linear temperature loading ( 1 0 T ≠ ), the thermal distributions of the in-plane and transverse shear stresses of a single-layer composite plate are shown in Figure 8. Similar to Figure 5, it can be observed that the EFSDTM_TN and LCW provide reliable solutions in describing the in-plane and transverse shear stresses of a single-layer composite plate under linear temperature loading.Furthermore, as shown in Figure 8, the accuracy of all of the theories considered in this study improved, relative to those obtained under uniform temperature loading.This tendency is attributed to the fact that linear temperature loading can cause bending behavior of the plates.For linear temperature loading (T 1 ̸ = 0), the thermal distributions of the in-plane and transverse shear stresses of a single-layer composite plate are shown in Figure 8. Similar to Figure 5, it can be observed that the EFSDTM_TN and LCW provide reliable solutions in describing the in-plane and transverse shear stresses of a single-layer composite plate under linear temperature loading.Furthermore, as shown in Figure 8, the accuracy of all of the theories considered in this study improved, relative to those obtained under uniform temperature loading.This tendency is attributed to the fact that linear temperature loading can cause bending behavior of the plates.The thermal distributions of the in-plane and transverse shear stresses of the sandwich plates under uniform temperature loading ( 0 0 T ≠ ) are shown in Figure 7.As shown in Figure 7a, the thermal solutions obtained by EFSDTM_TN and LCW are in good agreement with the exact solutions by precisely describing the severe noncontinuous distribution of the in-plane stress.In addition, Figure 7b indicates that EFSDTM_TN and LCW provide the best compromised thermal solutions for the local distribution of the transverse shear stress.Consequently, Figures 5-7 demonstrate that the transverse normal strain effect should be considered to accurately describe the thermal behavior of composite materials and sandwich structures subjected to uniform temperature loading.For linear temperature loading ( 1 0 T ≠ ), the thermal distributions of the in-plane and transverse shear stresses of a single-layer composite plate are shown in Figure 8. Similar to Figure 5, it can be observed that the EFSDTM_TN and LCW provide reliable solutions in describing the in-plane and transverse shear stresses of a single-layer composite plate under linear temperature loading.Furthermore, as shown in Figure 8, the accuracy of all of the theories considered in this study improved, relative to those obtained under uniform temperature loading.This tendency is attributed to the fact that linear temperature loading can cause bending behavior of the plates.[0 • /Core/0 • ] sandwich plates, respectively.From Figure 9, it can be concluded that the solutions obtained by EFSDTM_TN and LCW closely approximate the exact solutions, while the local solutions obtained from other theories are relatively inaccurate.Considering the thermal behavior of the sandwich plate, as shown in Figure 10, it is noteworthy that the EFSDTM_TN provides the best compromised solution for the inplane thermal stress by accurately capturing the noncontinuous local distribution.It can thus be concluded on the basis of Figures 3-10 that EFSDTM_TN and LCW can provide reliable thermomechanical solutions for laminated composites and sandwich plates because these theories reasonably consider the transverse normal strain effect.Although LCW provides the most accurate solution for some thermal problems, EFSDTM_TN has a prominent computational advantage due to its C0-based 5-DOF FE implementation, which can be highly compatible with commercial FE software.Therefore, it can be concluded that the FE implementation based on EFSDTM_TN is a useful approach in the thermomechanical analysis of laminated composites and sandwich plates.Core sandwich plates, respectively.From Figure 9, it can be concluded that the solutions obtained by EFSDTM_TN and LCW closely approximate the exact solutions, while the local solutions obtained from other theories are relatively inaccurate.Considering the thermal behavior of the sandwich plate, as shown in Figure 10, it is noteworthy that the EFSDTM_TN provides the best compromised solution for the in-plane thermal stress by accurately capturing the noncontinuous local distribution.It can thus be concluded on the basis of Figures 3-10 that EFSDTM_TN and LCW can provide reliable thermomechanical solutions for laminated composites and sandwich plates because these theories reasonably consider the transverse normal strain effect.Although LCW provides the most accurate solution for some thermal problems, EFSDTM_TN has a prominent computational advantage due to its C0-based 5-DOF FE implementation, which can be highly compatible with commercial FE software.Therefore, it can be concluded that the FE implementation based on EFSDTM_TN is a useful approach in the thermomechanical analysis of laminated composites and sandwich plates.Core sandwich plates, respectively.From Figure 9, it can be concluded that the solutions obtained by EFSDTM_TN and LCW closely approximate the exact solutions, while the local solutions obtained from other theories are relatively inaccurate.Considering the thermal behavior of the sandwich plate, as shown in Figure 10, it is noteworthy that the EFSDTM_TN provides the best compromised solution for the in-plane thermal stress by accurately capturing the noncontinuous local distribution.It can thus be concluded on the basis of Figures 3-10 that EFSDTM_TN and LCW can provide reliable thermomechanical solutions for laminated composites and sandwich plates because these theories reasonably consider the transverse normal strain effect.Although LCW provides the most accurate solution for some thermal problems, EFSDTM_TN has a prominent computational advantage due to its C0-based 5-DOF FE implementation, which can be highly compatible with commercial FE software.Therefore, it can be concluded that the FE implementation based on EFSDTM_TN is a useful approach in the thermomechanical analysis of laminated composites and sandwich plates.

Conclusions
In this study, an FE formulation based on EFSDTM_TN was developed and numerically validated for the reliable thermomechanical analysis of laminated composites and sandwich plates.The main features of the proposed FE model are summarized as follows: • MVT was employed in the proposed FE model to independently assume the displacement (FEST_TN) and transverse stress (EHOPT_TN) fields.The displacement and transverse stress fields were systematically interconnected in the MVT by establishing reasonable energy relationships.Based on the predefined relationships, the proposed FE model can not only embrace the explicit computational advantages of FSDT_TN, such as the C0-based 5-DOF FE implementation, but also ensure the solution accuracy of EHOPT_TN.• The transverse displacement field was enhanced by incorporating the components of external temperature loading to account for the contribution of the transverse normal strain effect efficiently.Consequently, the proposed FE model can provide reliable thermal solutions without introducing additional unknown variables.
In the proposed FE model, an 8-node serendipity element was employed to effectively derive higher-order derivatives while evaluating stress distributions.To demonstrate the numerical performance of the proposed FE model, several cases of thermalmechanical problems for laminated composites and sandwich structures were analyzed.The solutions obtained herein were then compared with those of conventional theories (FSDT, HSDT, and LCW), as well as 3-D exact solutions.From the numerical results, it can be concluded that the proposed FE model based on FESDTM_TN provides reliable thermomechanical solutions for laminated composites and sandwich plates.Consequently, it is expected that the proposed FE model can be applied to the thermomechanical analysis of laminated composites and sandwich structures with arbitrary geometries, loadings, and boundary conditions.[B] (m,b,s) membrane, bending, and transverse shear part of the strain matrix.
[K] e stiffness matrix for each element.
[K] (m, mb, b, s) membrane, membrane-bending coupling, bending, and transverse shear part of the element stiffness matrix.
[F] e M external force vector for each element derived from the mechanical loading.
[F] e T external force vector for each element derived from the thermal loading.
T mechanical loading applied in each element.
Thermo-Mechanical Problem 2.1.Mixed Variational Theorem Laminated composites and sandwich plates were considered as numerical models of thermomechanical problems.The geometric shapes and reference coordinates of the laminated plates are shown in Figure 1.Unless otherwise specified in the tensor notation, the Greek indices use values from set {1, 2}, whereas the Latin indices are assigned values from set {1, 2, 3}.L α and h represent the in-plane length and thickness of the laminated plates, respectively.x 3 denotes the transverse position which takes values within the range [−h/2, h/2].
. Unless otherwise specified in the tensor notation, the Greek indices use values from set {1, 2}, whereas the Latin indices are assigned values from set {1, 2, 3}.L α and h represent the in-plane length and thickness of the laminated plates, respectively.3 x denotes the transverse position which takes values within the range [

Figure 1 .
Figure 1.Geometric shape and reference coordinates of laminated plates.

Figure 1 .
Figure 1.Geometric shape and reference coordinates of laminated plates.

Figure 2 .
Figure 2. Comparison between analytical and FE solutions for stresses of [0 / 90 / 0 ] o o o laminated composite plates: (a) transverse shear stresses under mechanical loading; (b) in-plane stresses under uniform temperature loading; (c) transverse shear stresses under linear temperature loading.

Mathematics 2024, 12 , 963 14 of 21 Figure 3
Figure 3 shows the mechanical solutions of the in-plane displacements and stresses for cross-ply laminated composite plates.From Figure 3a, it can be seen that the EFSDTM_TN precisely describes the unsymmetrical zigzag distribution of the in-plane displacement in [0 / 90 / 0 / 90 ] o o o o laminated composite plates.As shown in Figure 3b, EFSDTM_TN can provide a reliable local solution for the in-plane stress by capturing its noncontinuous distribution.However, other theories are only useful for predicting the global behavior of in-plane displacement.

Mathematics 2024, 12 , 963 14 of 21 Figure 3
Figure 3 shows the mechanical solutions of the in-plane displacements and stresses for cross-ply laminated composite plates.From Figure 3a, it can be seen that the EFSDTM_TN precisely describes the unsymmetrical zigzag distribution of the in-plane displacement in [0 / 90 / 0 / 90 ] o o o o laminated composite plates.As shown in Figure 3b, EFSDTM_TN can provide a reliable local solution for the in-plane stress by capturing its noncontinuous distribution.However, other theories are only useful for predicting the global behavior of in-plane displacement.

Figure 5 .
Figure 5. Thermal solutions for a single-layer composite plate under uniform temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figure 5 .
Figure 5. Thermal solutions for a single-layer composite plate under uniform temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figure 5 .
Figure 5. Thermal solutions for a single-layer composite plate under uniform temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figure 8 .
Figure 8. Thermal solutions for a single-layer composite plate under linear temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figure 8 .
Figure 8. Thermal solutions for a single-layer composite plate under linear temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figure 8 .
Figure 8. Thermal solutions for a single-layer composite plate under linear temperature loading: (a) in-plane stresses; (b) transverse shear stresses.

Figures 9
Figures 9 and 10 illustrate the corresponding thermal behavior when the linear temperature loading ( 1 0 T ≠ ) is applied to [0 / 90 / 0 ] o o o

Figures 9
Figures 9 and 10 illustrate the corresponding thermal behavior when the linear temperature loading ( 1 0 T ≠ ) is applied to [0 / 90 / 0 ] o o o laminated composite and
αβ function of material properties to satisfy shear continuity conditions in the transverse stress field.properties to satisfy continuity conditions of transverse normal displacement.C N α , C M α in-plane correction factors derived by matching the force and moment resultants.

cN
,M(u * ,T 0 ,T 1 ) αβ coefficient of in-plane correction factors (C N α , C M α ).N i shape functions for the element in FE implementation.d e unknown displacements for each element.

F
(m, b, s) T o membrane, bending, and transverse shear part of the external force vector derived from the uniform temperature loading.F (m, b, s) T 1membrane, bending, and transverse shear part of the external force vector derived from the linear temperature loading.

Table 1 .
List of layup configurations for composite and sandwich plates.

Table 2 .
Convergence rate of central deflections for EFSDTM_TN.

Table 3 .
Kinematic unknown variables and total DOFs of the FE models (32 × 32 mesh density).

Table 3 .
Kinematic unknown variables and total DOFs of the FE models (32 × 32 mesh density).