Study on the effects of initial stress and imperfect interface on static and dynamic problems in thermoelastic laminated plates

1 Introduction

With the rapid development of modern engineering technologies, composite materials have been widely used in aerospace, civil engineering, and automotive manufacturing fields due to their unique physical and chemical properties. Among them, composite laminates, as a typical structural form, have long been a research focus in the fields of materials science and mechanics regarding their thermoelastic problems. Existing studies have shown that the coupling between temperature and strain fields not only regulates the propagation speed of waves but also affects their waveform, leading to the attenuation of elastic waves [1]. Therefore, in complex and varied engineering environments, evaluating the mechanical performance and safety of laminate structures under thermoelastic coupling conditions, as well as constructing an accurate and efficient analytical model for thermoelastic materials, is crucial for ensuring the stability and reliability of composite laminate structures.

Due to limitations imposed by mechanical loads, thermal fields, and manufacturing techniques, initial stresses in layered structures are inevitable in practical engineering applications. Consequently, the issue of wave propagation in layered structures under initial stresses has gained increasing attention. Abd-Alla et al. [2] employed a bisection method to investigate the influence of initial stresses on the dispersion relation of transverse waves in piezoelectric layered cylinders. Cilli and Ozturk [3] investigated the dispersion of torsional waves in layered circular cylinders with initial stresses. Based on the principles of incremental deformation mechanics, Zhang and Xue [4] utilized the Legendre polynomial series method to examine guided waves in layered hollow cylinders in both radial and axial directions. Du et al. [5] conducted a theoretical analysis of SH wave propagation in circular cylindrical layered piezoelectric structures with initial stresses. Zhao et al. [6] studied the SH waves propagating in layered cylinders under initial stresses. Based on the 3D linearized theory, Akbarov et al. [7] studied the impact of initial strains on the dispersion of torsional waves in hollow sandwich cylinders. Yu et al. [8] adopted a double-orthogonal polynomial method to analyze the effect of axial initial stresses on wave characteristics in layered sectorial cylindrical structures, calculating the displacement distributions and dispersion curves for various sectorial cylindrical structures.

For plate structures, many literature sources have studied the effect of initial stress on the propagation of bulk or surface waves in layered media. Liu et al. [9] investigated the propagation of Love waves in layered piezoelectric structures with initial stress, discussing the effects of uniform initial stress on phase velocity, stress field, and electromechanical coupling factors. Qian et al. [10] and Su et al. [11] studied the influence of non-uniform initial stress on Love waves in different layered structures, respectively, revealing that initial stress has significant effects on dispersion relations, phase velocity, group velocity, and other parameters. Singh [12] studied wave propagation in prestressed piezoelectric half-spaces, obtaining reflection coefficients for qP and qSV waves. Kayestha et al. [13] considered the dispersive behavior of time-harmonic waves propagating along a principal direction in perfectly bonded bi-material laminated plates. Du et al. investigated the propagation of Love waves [14] and SH waves [15] in magneto-electro-elastic laminated structures with uniform initial stress. Zhang et al. [16] investigated the influence of inhomogeneous initial stress on the propagation behavior of Love waves in multi-field layered media. Yang [17] studied the effect of biasing field stress on the free vibration in an electroelastic body.

In the previous research on layered structures, it is assumed that the interlayer bonding is perfect, meaning that the field variables at the interfaces satisfy continuity conditions. However, in practical engineering, manufacturing techniques or the aging of adhesive layers may degrade the interfacial bonding strength, leading to interfacial slipping or debonding. The influence of such bonding defects on structural behavior has been extensively studied [18,19,20], with most of these studies based on various simplified plate and shell theories. Recently, Serpilli et al. utilized the asymptotic expansion techniques to model bonded interfaces in linear multiphysics composites, proposing the soft interface models, hard interface models, and rigid interface models [21], and applied these techniques to the piezoelectric and thermo-mechanical coupling frameworks [22]. Serpilli et al. [23] also introduced several novel interface conditions between layers of three-dimensional composite structures within the framework of coupled thermoelasticity. In subsequent work, Serpilli et al. [24] studied the mechanical behavior of composites made of piezoelectric hollow spheres under both full and incomplete contact conditions using the asymptotic method and the transfer matrix homogenization method, and validated the correctness and effectiveness through numerical examples. Nevertheless, Chen and Lee [25] revealed that directly applying traditional higher-order theories to layered structures with imperfect interfaces may result in significant errors. Wang and Pan [26] employed a special form of transfer matrix to provide exact solutions for simply supported multi-layered piezoelectric plates with thermal effects and imperfect interfaces. Chen et al. [27] studied the bending problem of multiferroic rectangular plates with imperfect interfaces through a 3D exact theory. A generalized spring-layer model was proposed to characterize imperfect interfaces and establish a state-space formulation incorporating an interfacial transfer matrix, which can be conveniently applied to piezoelectric, piezomagnetic, and elastic laminated plates.

From the brief review above, it is evident that extensive research was accomplished on wave propagation in layered cylindrical structures with initial stress, but relatively few studies have analyzed the vibration and static bending of layered thermoelastic structures under initial stress, especially concerning the impact of interfacial defects and initial stress on thermoelastic coupling media. Therefore, this study extends upon our previous work [28] by investigating the dynamic and static problems of thermoelastic laminated plates under cylindrical bending, comprehensively considering the effects of initial stresses and interface defects in both the elastic and thermal fields. The imperfect bonding is modeled via a spring-layer model, and two types of conditions in the thermal field are adopted similar to those in previous studies [25–28]. By introducing the interfacial transfer matrix into the state-space analysis, the effects of interfacial defects can be easily taken into account. Numerical examples demonstrate that different interfacial defect conditions and initial stresses have significant effects on the natural frequencies and modes of the structure, which deserve attention. Furthermore, the numerical analysis of static bending also reveals the impacts of the degree, type of interfacial defects, and initial stress, which are crucial parameters to consider.

2 Basic equations

Consider an N-layered cross-ply thermoelastic laminated plate under the assumption of cylindrical bending deformation, as depicted in Figure 1. In this case, there are only two non-vanishing displacements u and w, along x and z directions, respectively. The thermoelastic body is bearing initial stresses. The temperature difference on the top layer of the plate is T ₁, while the temperature difference on the bottom layer is T ₀.

Figure 1

Multi-layered laminated plate.

The constitutive and heat conduction relations, the mechanical equilibrium, and heat energy equations for a thermoelastic laminated plate under cylindrical bending are as follows:

where T is the temperature difference relative to a referenced temperature at the stress-free state; σ i and τ i j the normal and shear stresses; u and w the displacements along x and z directions, respectively; P i ( i = x , z ) are the heat flux; c i j , β i , and k ii , respectively, the elastic, thermal elastic modulus, and heat conduction coefficients; σ i 0 ( i = x , z ) the initial stresses.

Considering the cylindrical bending problem, from Equations (1) and (2), it is straightforward to derive the following state equation as refs. [29,30]:

where the operator matrix U is given in Equation (A1) of Appendix, F = [ u σ z T τ x z w P z ] T are the six state variables, in terms of which the two induced variables, σ x and P x , can be expressed as Equations (A2) and (A3) in Appendix.

3 Exact solution

Consider the particular case of a thermoelastic laminate in cylindrical bending with simple support, for which an exact solution can be obtained.

For a thermoelastic laminate, the simply supported conditions can be described:

Compared with an elastic plate [31], the thermal condition is as follows:

To satisfy the aforementioned boundary conditions, one can assume the following:

where ζ = z / h and ξ = x / l are the dimensionless coordinate variables, n is the wave number, ω is the circular frequency, and parameters with superscript (1) are the material constants of the first (bottom) layer, as depicted in Figure 1.

Substituting Equation (6) into Equation (3) yields [28]

where F = [ u ¯ σ ¯ z T ¯ τ ¯ x z w ¯ P ¯ z ] T is the dimensionless state vector, and the expression of the coefficient matrix U e is given in Appendix (A4). The solution to Equation (7) is

where ζ 0 = 0 , ζ i = z i / h = Σ j = 1 i h j / h , and h i is the thickness of the ith layer. Let ζ = ζ k in Equation (8) yield

If one replaces sin n π ξ by cos n π ξ , and vice versa, in Equation (6), another exact solution could be similarly derived. In this case, the boundary condition at the end is known as the guided condition (or rigidly slipping or rigidly smooth contact) [32], i.e.,

Accordingly, u could be proposed as sinusoidal or cosinoidal.

4 Interfacial transfer matrix

For an arbitrary interface, such as the one between the kth and (k + 1)th layer ( z = z k ), imperfect mechanical conditions follow the spring-layer model [25,26,27,28]:

In [26,33], the authors proposed two imperfect interfaces in heat conduction: the weak thermal conducting condition,

and the high thermal conducting condition,

where R x x ( k ) , R z z ( k ) , R w ( k ) , and R h ( k ) represent the elastic stiffness coefficients, the weak and high interface coefficients for the thermal field, respectively. Further discussions on the mechanical spring-layer model are reported in previous studies [34,35]. This model can be analogously applied to electric, magnetic, or thermal fields. By introducing the imperfect interface model, the influence of imperfect interfaces in the structure can be readily taken into consideration.

Using Equations (6), (11) and (12) or (13) can be expressed in the following matrix form:

where P k represents the interfacial transfer matrix for weak and high thermal conducting problems, respectively,

with R ¯ x x ( k ) = c 55 ( 1 ) R x x ( k ) h , R ¯ z z ( k ) = c 55 ( 1 ) R z z ( k ) h , R ¯ w ( k ) = k 33 ( 1 ) R w ( k ) h , R ¯ h ( k ) = − R h ( k ) ( n π s ) 2 k 33 ( 1 ) h .

The relationship [36] between the state vectors of the top and bottom surfaces of the laminated plate can be readily obtained from Equations (9) and (14):

where F 1 ( N ) and F 0 ( 1 ) are the state vectors at the top and bottom surfaces, respectively, and H represents the global transfer matrix for a laminate with imperfect interface, and P 0 = I . In the case of perfectly bonded laminated plates, all P j become units, yielding H = ∏ j = N 1 M j .

One usually has the following mechanical boundary conditions:

However, the temperature conditions are as follows:

Substitution of Equations (18) and (19) into Equation (17) determines the unknown variables on the bottom layer:

where H i j are the elements of the global transfer matrix H . Once all the unknowns on the bottom layer are solved, the variables at an arbitrary point can be derived using Equations (9) and (14).

5 Numerical examples

Based on the previous derivation, this section will provide several numerical examples of composite plates with imperfect interfaces and initial stresses, including natural frequencies, vibration modes, and static bending distributions.

First, to verify the correctness of the formula derivation and the program, the analytical results from the previous sections are compared with the numerical results. A single-layer rectangular thin plate with simply supported boundaries at both ends is selected. The geometric parameters are thickness h = 2 mm , width l = 20 mm , length = 184.4 mm, density ρ = 7.85 × 10³ kg/m³, elastic modulus E = 200 GPa, and Poisson’s ratio 0.3. Both the initial stress and the weak interface coefficient are set to zero in this example. The natural frequencies of the simply supported plate for free vibration are calculated, and the analytical results, finite element results, and their errors are shown in Table 1.

As seen from Table 1, the analytical results obtained through the aforementioned theoretical derivation show small errors compared to the finite element simulation results, with the maximum error not exceeding 10%, thereby verifying the reliability of the derivation and program. The errors arise from the fact that the finite element model can only represent a finite-length plate, rather than an infinite-length plate, thus leading to some discrepancies with the analytical solution for plane strain problems. Therefore, using this model, one can effectively analyze the static and dynamic characteristics of thermoelastic laminated plates, as well as the influence of initial stresses and interface defects on various field variables.

For all the following numerical examples, the material parameters considered are listed in Table 2 [37], where material I is Si₃N₄ and material II is Cobalt. It is assumed that the thickness of each layer in the two-layer plate is the same, and the dimensionless aspect ratio s = h / l and dimensionless frequency Ω 2 = ρ ( 1 ) ω 2 h 2 / c 55 ( 1 ) are given. The imperfect interface coefficients are R ¯ x x ( k ) = R ¯ z z ( k ) = R ¯ w ( k ) = R ¯ h ( k ) = R . Hereinafter, the mechanical boundary conditions are assumed to be traction-free at the top and bottom surfaces. Assume that the temperature difference on the bottom surface is T 0 = c 55 ( 1 ) T ¯ 0 / β 3 ( 1 ) , and the temperature difference on the top surface is T 1 = c 55 ( 1 ) T ¯ 1 / β 3 ( 1 ) .

5.2 Modal analysis

This section investigates the effect of initial stress and imperfect interfaces on vibration modes. Figures 2 and 3 show the modal shapes of various field vectors in the z-direction under conditions of weak thermal conductivity. A two-layer structure is considered, with material properties listed in Table 2 and other parameters being s = 0.1, R = 0.05, 0.1, 0.2, σ x 0 = σ z 0 = 10 8 α , and α = 1, 0, −1.

Figure 2

Modal analysis of the laminated plate under weak thermal conducting condition ([I/II], n = 3, R = 0.05).

Figure 3

Modal analysis of the laminated plate under weak thermal conducting condition ([I/II], n = 1, α = 0 ).

Figure 2 presents the modal shapes under different initial stress conditions, studying how various field variables change along the z-direction as the initial stress varies. In this example, three cases with initial stress coefficients α of 1, 0, and −1 are selected. The blue solid, red dashed, and black dash-dotted lines in the figure correspond to the results of α = 1 , α = 0 , and α = − 1 , respectively. These three curves represent the modes with dimensionless frequencies of 0.29077502, 0.28988029, and 0.28898282 (wave number n = 3). As seen from Figure 2, the modal curves for u , σ z , τ x z , and w overlap, but this is not the case for T and P z . For α = 0 and 1, the modal curves for T and P z overlap; however, for α = − 1 , T , and P z have opposite values compared to the other two cases.

To analyze the influence of interface defects on various field variables, Figure 3 shows the changes in field variables along the z-direction under different imperfect interface coefficients. Three cases with the imperfect interface coefficients of R = 0.05, 0.1, and 0.2 are selected in Figure 3. The blue solid, red dashed, and black dash-dotted lines correspond to the results for R = 0.05, 0.1, and 0.2, respectively. These three curves represent the modes with dimensionless frequencies of 0.03780132, 0.03774749, and 0.03764085 (n = 1). As can be seen from Figure 3, different imperfect interface coefficients have a significant impact on T and P z ; especially, the curve for R = 0.1 has opposite signs compared to the curves for R = 0.05 and R = 0.2. The modal curves for other field variables overlap.

5.3 Static bending

This section investigates the effects of initial stress, imperfect interfaces, and thermal boundary conditions on the field vectors for static bending. A two-layer structure is considered, with material properties detailed in Table 2 and other parameters being s = 0.1, n = 1, T ¯ 1 = 2 , T ¯ 0 = 1 , σ x 0 = σ z 0 = 10 8 α , and α = 1, 0, −1.

To study the effect of different initial stresses on field variables under static bending conditions, Figure 4 presents the variation of field variables along the z-direction as the initial stress changes under conditions of weak thermal conductivity. Figure 4 presents three cases with initial stress coefficients α of 1, 0, and −1, where the blue solid, red dashed, and black dash-dotted lines correspond to the results of α = 1 , 0, and −1, respectively. As seen from Figure 4, only the distribution plots of w show slight differences among three different initial stresses, while those for other variables almost overlap. Additionally, the temperature T has discontinuities at the interface ( ζ = 0.5 ), while other quantities have sharp corners at the interface, indicating that it is not differentiable at these points.

Figure 4

Static bending of the laminated plate under weak thermal conducting condition ([I/II], n = 1, R = 0.05).

To further analyze the influence of the imperfect interface coefficients on the static bending, Figure 5 presents the distribution of field variables along the z-direction, for the weak thermal conductivity condition. In Figure 5, the blue solid, red dashed, and black dash-dotted lines correspond to the results for three imperfect interface coefficients of R = 0.05, 0.1, and 0.2, respectively. As seen from Figure 5, there are relatively small differences in u , σ z , and τ x z , but significant differences in T , w , and P z . As the imperfection of the interface worsens, the change of T and P z increases, while the amplitude of w decreases. Except for u , σ z , and P z , the remaining variables exhibit more pronounced jumps at the interface ( ζ = 0.5 ).

Figure 5

Static bending of the laminated plate under weak thermal conducting condition ([I/II], n = 1, α = 0 ).

Figure 6 compares the trends of field variables along the z-direction for conditions of high and weak thermal conductivity. The blue solid line represents the weak thermal conductivity condition, while the red dashed line represents the high thermal conductivity condition. As can be seen from Figure 6, differences in temperature T are observed under the weak thermal conductivity condition, while differences in heat flux P z are observed under the high thermal conductivity condition, which is consistent with the assumption of Equations (12) and (13). However, the other quantities show little difference under the two thermal conditions. For both thermal conditions, there are discontinuities or abrupt changes of stresses τ x z , displacement w , and temperature T at the interface ( ζ = 0.5 ), indicating the possibility of delamination at the laminate interface.

Figure 6

Static bending of the laminated plate under weak/high thermal conducting conditions ([I/II], n = 3, R = 0.05, α = 1 ).

6 Conclusion

Based on the generalized thermoelastic theory, this study analyzes the influence of initial stresses and interfacial imperfections in thermoelastic laminated plates, yielding analytical solutions for two-dimensional thermoelastic rectangular laminates. Beyond the spring-layer model for the elastic field, an analogous linear relationship is introduced for temperature differences and heat flux in the thermal field. By incorporating an interfacial transfer matrix and interlayer transfer matrix that encapsulate imperfect elastic and thermal interface conditions, various imperfect interface issues in multi-layer thermoelastic materials can be discussed within a state-space framework, encompassing two heat conduction conditions: weak and high thermal conductivity. Further parametric analysis leads to the following conclusions:

1. As evident from Table 3, the dimensionless frequency increases with the augmentation of the initial stress coefficient. Specifically, initial tensile stress enhances the dimensionless frequency, whereas initial compressive stress reduces it.

2. Regardless of whether the initial stress in the structure is zero or non-zero, an increase in the imperfect interface coefficient leads to a decrease in the dimensionless frequency. The presence of weak interfaces degrades stiffness, subsequently lowering the frequency.

3. Initial stresses and weak interfaces have minimal impacts on various modes of the elastic field but significantly influence temperature differences and heat flux.

4. Initial stresses slightly affect the distribution of displacement w under static bending. Weak interfaces, on the other hand, have a profound influence on the distributions of temperature differences T , displacements w , and heat flux P z for static bending. Slight differences in displacement w , heat flux P z , and temperature differences T are observed under different thermal conditions. Discontinuities or abrupt changes of stresses, displacements, and temperatures at the interface suggest the possibility of interface debonding.

These findings may contribute to elucidating the behavior of practical thermoelastic laminated structures, particularly in health monitoring for such structures. Our solutions can serve as benchmarks for various plate theories and approximation methods, as well as for calibrating the outcomes of other numerical methods.