Stress-Function Variational Method for Accurate Free-Edge Interfacial Stress Analysis of Adhesively Bonded Single-Lap Joints and Single-Sided Joints

Large free-edge interfacial stresses induced in adhesively bonded joints (ABJs) are responsible for the commonly observed debonding failure in ABJs. Accurate and efficient stress analysis of ABJs is important to the design, structural optimization, and failure analysis of ABJs subjected to external mechanical and thermomechanical loads. This paper generalizes the high-efficiency semi-analytic stress-function variational methods developed by the authors for accurate free-edge interfacial stress analysis of ABJs of various geometrical configurations. Numerical results of the interfacial stresses of two types of common ABJs, i.e., adhesively bonded single-lap joints and adhesively single-sided joints, are demonstrated by using the present method, which are further validated by finite element analysis (FEA). The numerical procedure formulated in this study indicates that the present semi-analytic stress-function variational method can be conveniently implemented for accurate free-edge interfacial stress analysis of various type of ABJs by only slightly modifying the force boundary conditions. This method is applicable for strength analysis and structural design of broad ABJs made of multi-materials such as composite laminates, smart materials, etc.


Introduction
High-performance polymeric adhesives are commonly utilized for connecting separate parts in engineering practices to realize their structural integrity and functions such as load transfer, stiffness enhancement, surface repairing, etc., which has resulted in various adhesively bonded joints (ABJs), as illustrated in Figure 1. Compared to their counterparts of mechanically-fastened bolted, riveted, and welded joints, ABJs carry unique structural and mechanical advantages such as simplified structural design and fabrication, reduced joining space and weight, enhanced fatigue durability, extended crack growth tolerance, suppression of noises, and so on [1][2][3]. So far, advanced adhesive joining techniques have been integrated into manufacturing of various load-carrying structures in modern aircrafts, marine and ground vehicles, etc. [4][5][6]. For example, adhesively bonded metallic joints have been successfully structured in commercial aircrafts with the advent of Airbus A300 and modern Boeing aircrafts (e.g., Boeing 737) [7][8][9]. Besides extensive deployment of ABJs in structural applications, adhesive joining technology has also been broadly utilized in microelectronics packaging since 1970s. In fabrication of electronic devices, it is the common practice to join multiple heterogeneous materials together by adhesives or solders. Temperature increase due to heat generation in electronic devices during service commonly induces high edge interfacial stresses between bonded materials, resulting in structural failure and function degradation. Substantial theoretical and experimental studies have been conducted in the last four decades for accurate and efficient determination of the interfacial thermomechanical stresses in bonded thermostats of microelectronic devices [10][11][12][13][14][15][16][17][18][19], in In view of solid mechanics, exact stress analysis of ABJs is mathematically challenging even in the simple cases of linearly elastic adherends and adhesives due to the complexity of solving a set of coupled partial differential equations (PDEs) with multiple boundary conditions (BCs). Historically, several classic joint models have been proposed for determining the stress distribution along the bonding lines and related debonding failure. Nearly a century ago, Timoshenko [20] was the first to realize the stress concentration near free-edges of a bimaterial thermostat. Volkersen [21] and Goland and Reissner [22] were the pioneers who obtained the first linear elasticity solutions to the interfacial stresses in single-lap ABJs, and their explicit interfacial stress solutions are widely regarded as the classic stress solutions of ABJs. These stress solutions have been extensively adopted in various textbooks and in almost all research papers on stress analysis of ABJs late on. In these ABJ models, a simple assumption is made such that the adhesive layers are treated as one-dimensional (1D) pure shear springs (i.e., the shear-lagging adhesive model) and they do not obey the generalized Hooke's law as only the shear modulus is involved. Such an oversimplified shear-lagging assumption has been widely adopted by most late investigators in finding the interfacial stresses in various types of ABJs. Yet, obvious limitations have been identified in these models and related stress results: The peak shear stresses occur at the free adherend ends that obviously violate the shear-free condition at the free-ends, and noticeable stress variations across the adhesive layer near free edges are ignored, among others [2,[23][24][25][26].
With aid of the shear-lagging assumption of the adhesive layers in ABJs by Volkersen [21] and Goland and Reissner [22], quite a few elegant analytic solutions of interfacial stresses have been obtained for various ABJs with different extents of deliberation in the past four decades. To mention a few, Erdogan and Ratwani [27] formulated a generalized 1D ABJ model for prediction of the interfacial stresses in adhesively stepped joints under uniaxial tension, which leads to closed-form solutions of the normal and shear interfacial stresses of the joints. In the limiting case of an infinite number of steps, this joint model can be used to predict the interfacial stress variations in adhesively bonded scarf joints under uniaxial tension. In this ABJ model, the shear-lagging assumption of the adhesive layers plays a critical role in coupling the governing equations of the individual adherends. Based on the same assumption, Delale et al. [28] further formulated the ABJ models to universally determine the stress field in adhesively bonded single-lap joints and single-sided strap joints, in which the adherend layers are treated as elastic plates under In view of solid mechanics, exact stress analysis of ABJs is mathematically challenging even in the simple cases of linearly elastic adherends and adhesives due to the complexity of solving a set of coupled partial differential equations (PDEs) with multiple boundary conditions (BCs). Historically, several classic joint models have been proposed for determining the stress distribution along the bonding lines and related debonding failure. Nearly a century ago, Timoshenko [20] was the first to realize the stress concentration near free-edges of a bimaterial thermostat. Volkersen [21] and Goland and Reissner [22] were the pioneers who obtained the first linear elasticity solutions to the interfacial stresses in single-lap ABJs, and their explicit interfacial stress solutions are widely regarded as the classic stress solutions of ABJs. These stress solutions have been extensively adopted in various textbooks and in almost all research papers on stress analysis of ABJs late on. In these ABJ models, a simple assumption is made such that the adhesive layers are treated as one-dimensional (1D) pure shear springs (i.e., the shear-lagging adhesive model) and they do not obey the generalized Hooke's law as only the shear modulus is involved. Such an oversimplified shear-lagging assumption has been widely adopted by most late investigators in finding the interfacial stresses in various types of ABJs. Yet, obvious limitations have been identified in these models and related stress results: The peak shear stresses occur at the free adherend ends that obviously violate the shear-free condition at the free-ends, and noticeable stress variations across the adhesive layer near free edges are ignored, among others [2,[23][24][25][26].
With aid of the shear-lagging assumption of the adhesive layers in ABJs by Volkersen [21] and Goland and Reissner [22], quite a few elegant analytic solutions of interfacial stresses have been obtained for various ABJs with different extents of deliberation in the past four decades. To mention a few, Erdogan and Ratwani [27] formulated a generalized 1D ABJ model for prediction of the interfacial stresses in adhesively stepped joints under uniaxial tension, which leads to closed-form solutions of the normal and shear interfacial stresses of the joints. In the limiting case of an infinite number of steps, this joint model can be used to predict the interfacial stress variations in adhesively bonded scarf joints under uniaxial tension. In this ABJ model, the shear-lagging assumption of the adhesive layers plays a critical role in coupling the governing equations of the individual adherends. Based on the same assumption, Delale et al. [28] further formulated the ABJ models to universally determine the stress field in adhesively bonded single-lap joints and single-sided strap joints, in which the adherend layers are treated as elastic plates under cylindrical bending. Yet, the shear stresses predicted by this model do not satisfy the shear-free condition at the adherend ends [2,3]. Refined finite element analysis (FEA) indicated that the interfacial stresses predicted by this model are overshot in a large region from the adherend ends [29]. Chen and Cheng [30] proposed an ABJ model for stress analysis of adhesively bonded single-lap joints by assuming that the axial stresses in adherends vary linearly across the thickness (i.e., Euler-Bernoulli beam), and the shear stress across the adhesive layer is assumed to be constant. In this model, the entire stress field in the ABJ can be expressed in terms of two unknown normal stress functions via triggering the stress equilibrium equations within the framework of two-dimensional (2D) elasticity. These two unknown stress functions can be further determined via solving a set of two coupled fourth-order ordinary differential equations (ODEs) according to the principle of minimum complementary strain energy. The stress field gained in this ABJ model is able to satisfy all the traction BCs, and the predicted location of the peak interfacial shear stress appears at a distance of~20% of the adherend thickness from the adherend ends as validated quantitatively by refined FEA [29]. Yet, due to oversimplification of the adhesive layer, this ABJ model yielded a physically questionable zero normal stress in the adhesive layer along the bonding line. In addition, by using a simple shear-lagging model of the adhesive layer, Her [31] obtained the closed-form solutions to the axial force in the adherends and the shear-force in the adhesive layer of adhesively bonded single/double-lap joints, which were largely validated numerically by FEA. Yet, the static equilibrium in terms of bending moments of the joint and the shear-free conditions at the adherend ends were not satisfied. Thus, Her's solutions can only be restricted in the case of ABJs made of very thin adherends, in which the thickness effect is treated as the higher order terms and therefore the bending effect is ignored. In addition, Tsai et al. [13] extended the classic ABJ models formulated by Volkersen [21] and Goland and Reissner [22] via adopting the linear shear deformation across the adhesive layer. This model is able to recover the classic Volkersen's and Goland and Reissner's models in the limiting cases while it does not satisfy the shear-free conditions at adherend ends. By modeling the adhesive layer as two distributed linearly elastic shear and tension springs, Lee and Kim [32] derived the closed-form solutions to the axial force in adherends and the shear force in the adhesive layer of adhesively bonded single-lap joints (ABSLJs), which were validated by their detailed FEA, except for the fact that the fundamental shearfree conditions at the adherend ends were not satisfied. Radice and Vinson [33] formulated a higher-order ABJ model, in which Airy stress potential for a 2D elastic ABJ body is expressed as the sum of a series of power functions with respect to the thickness coordinate and is consequently determined via solving the resulting Cauchy-Euler equations in favor of Rayleigh-Ritz minimization of the potential energy of the entire ABJ.
So far, developing efficient and robust ABJ models for accurate stress analysis is still an active research topic that attracts researchers all over the world. A number of recent analytic solutions for strength and fracture analyses of traditional ABJs and layered materials have been formulated such as for composite and heterogeneous adherends, asymmetric joints, functionally gradient adhesive layers, etc. [3,. It needs to also be mentioned that in the theoretical approach, the stress singularity exponent of the free-edge stresses of ABJs depends upon the material properties of bonded materials, i.e., two Dundurs' parameters and the edge angle [59][60][61][62], and differs from that of interfacial cracks. More recent studies in design and structural reinforcement of ABJs also include development of smart ABJs integrated with piezoelectric layers [63][64][65][66] and toughening and damage self-healing nanofiber interlayers [67][68][69][70][71][72][73][74][75][76][77] for damage sensing and debonding suppression, in which the nonwoven continuous monolithic and core-shell nanofiber interlayers are produced by means of the low-cost top-down electrospinning technique [78][79][80][81][82][83][84][85][86]. Moreover, several layerwise joint models have also been formulated for improving the stress analysis of ABJs. For instance, Hadj-Ahmed et al. [87] proposed a layerwise ABJ model for a multi-layered ABJ that was modeled as a stack of Reissner plates to be coupled through the interlaminar normal and shear stresses. A set of governing ODEs was obtained via minimization of the total strain energy of the ABJ. Besides, Diaz et al. [88] further formulated an improved layerwise ABJ model, in which the ABJ was modeled as a stack of Reissner-Mindlin plates. As a result, a set of eight governing ODEs was extracted via evoking the constitutive laws and solved to satisfy the traction BCs. This ABJ model can be well validated by FEA for free-edge interfacial stress prediction. Moreover, Yousefsani and Tahani [36][37][38] proposed another version of the layerwise ABJ model. In this model, displacements of artificially divided sub-layers of the ABJ were treated as the field variables, and a set of governing ODEs was obtained via minimization of the total potential energy of the ABJ. For the purpose of accurate interfacial stress prediction, 18 artificial sub-layers were used in their numerical examples. Such layerwise ABJ models were further extended for stress analysis of smart joints integrated with piezoelectric patches [63,64]. Detailed literature surveys on the historical progress in mechanics of ABJs can be found in several recent review papers [6,[89][90][91][92][93][94][95][96] and references therein.
On the other hand, to effectively approach the stress conditions in ABJs, in particular the traction-free conditions at the free-edges of ABJs, Chang [97][98][99][100] expressed the interfacial peeling and shear stresses on the bonding lines in terms of the sums of an infinite series of sine or cosine functions, respectively, with their coefficients determined via minimization of the strain energy of the entire ABJs. During the process, the axial stress in each elastic adherend and adhesive layer was assumed to linearly vary across the thickness of the corresponding layer as that of classic Euler-Bernoulli beams, and the related transverse normal stresses and shear stresses were determined by evoking the 2D stress equilibrium equations. To simplify the process, Chang adopted the deformation (deflection) compatibility of bonded dissimilar adherends in bending within the framework of Euler-Bernoulli beam theory of composite beams. The advantages of Chang's approaches are that all the interfacial stress solutions can be expressed as the sums of infinite trigonometrical series, which can be further added up into elegant closed-form expressions. Yet, refined FEA indicated that Chang's approach noticeably underestimates the shear stress variation near the free edges of ABJs, due mainly to the harsh treatment of the deformation compatibility [23][24][25][26]29].
To overcome the above theoretical obstacles in stress analysis of ABJs within the classic Euler-Bernoulli beam theory, Wu and co-workers [2,3,[23][24][25][26]29,101,102] formulated an efficient stress-function variational method for accurate interfacial stress analysis of a variety of ABJs including bonded joints and adhesively bonded monolithic and composite joints. Different from other methods available in the literature, Wu's approaches introduce two unknown interfacial normal (peeling) and shear stress functions at each interface, and the axial stresses in the adherends and adhesive layers are assumed to linearly vary across the thickness as that of the classic Euler-Bernoulli beams. By evoking the 2D stress equilibrium equations, the rest planar stress components in the ABJs are expressed exactly in terms of the unknown interfacial stress functions at the upper and lower interfaces [2,23]. It was shown that such treatment guarantees all the stress components to be consistent across the bonding lines [2], which endorses that this method can be extended to determine stresses in multi-layered ABJ systems such as composite joints [101]. Finally, these unknown interfacial stress functions are determined via solving a set of coupled ODEs, which result from minimization of the complimentary strain energy of the entire joints. In the simple case of bonded joints made of two elastic adherends, a set of two coupled ODEs with respect to two interfacial stress functions can be obtained [23,25,26], which can be further reduced into one governing ODE via introducing the deformation compatibility of the two adherends within the classic Euler-Bernoulli beam theory [24]. In the case of ABJs made of two adherends adhesively bonded through an adhesive layer, a set of four coupled ODEs can be obtained with respect to two pairs of interfacial peeling and shear stress functions at two interfaces [2]. The interfacial peeling and shear stresses of the ABJs determined by this method can exactly satisfy the traction-free conditions at the free edges of the adherends. In addition, detailed FEA indicates the high accuracy of this semi-analytic stress-function variational method for stress analysis of ABJs [2,23].
With the above literature review, it can be concluded that a number of ABJ models with varying extents of deliberation have been formulated in the literature for stress and strength analysis and structural design of ABJs since the pioneering studies by Volkersen [21] and Goland and Reissner [22]. Yet, many such ABJ models overlooked or ignored the particular restrictions of using various assumptions, in particular the shear-lagging assumption, to conduct the stress analysis of ABJs. Thus, additional studies are still needed to elucidate such fundamental and important analysis. Specifically, in many existing ABJ models, (1) conditions of the shear-lagging assumption are typically broken due to the large shear and normal stresses at free edges; (2) shear and normal stresses have noticeable changes across the thin adhesive layer near the free edge, thus many ABJ models based on the assumption that stress variation across the adhesive layer is negligible are not self-consistent. In fact, no interfacial stresses at both the upper and lower surfaces of the adhesive layers are determined in those ABJ models to support and verify the assumptions adopted in those ABJ models; (3) many ABJ models predict the interfacial shear stresses that do not satisfy the simple shear-free conditions at the free edges. Nevertheless, among a few others, the recent ABJ models based on the semi-analytic stress-function variational method formulated by the authors showed the advantages of fundamentally resolving the above three issues and their accuracies were validated via detailed FEA. Thus, in this paper, we further generalize this effective and high-efficiency semi-analytic stressfunction variational method for determining the interfacial stresses in adhesively bonded single-lap joints (ABSLJs) and adhesively single-sided joints (ASSJs) to show its efficiency, accuracy, robustness, and universality for a broad range of ABJs. Detailed derivations of this generalization are given. Numerical examples and scaling analysis of the interfacial stresses of these two ABJs are conducted and compared. Discussions and conclusions of the present study are made in consequence.

Static Equilibrium Equations of General ABJs
As shown in Figure 1, a typical ABJ consists structurally of three homogeneous, isotropically linearly elastic segments, i.e., one three-layered joint segment and two single-layered adherend segments, though more structurally complicated ABJs can also be generated via attaching additional adhesive and adherend layers. The ABJ deformations are assumed small, and nonlinear geometrical effects are not considered (such as large deflections). In general, the two slender adherend segments can be treated simply as two Euler-Bernoulli beams, while accurate stress analysis is challenging in the three-layered joint segment due to its multiple materials and boundaries subjected to external loads. This three-layered segment is the focus in ABJ modeling that has attracted a large number of investigators in the past several decades. Herein, the study starts with the accurate stress analysis of this three-layered joint segment subjected to external mechanical and thermomechanical loads. In a general case of ABJs, due to loss of the lateral symmetry, deformations of this three-layered segment in an ABJ are a combination of axial elongation and lateral deflection. The upper/lower adherend and adhesive layers of the joint are assumed to be slender, with the same width b and the thicknesses of h 1 , h 2 , and h 0 , respectively, and their axial stresses can be approximated to follow the classic Euler-Bernoulli beam theory while the shear and lateral normal stresses are determined according to the static equilibrium equations in 2D elasticity. Free-body diagrams (FBDs) of the representative segments of the upper/lower adherend and adhesive layers are shown in Figure 2a-c, in which the stress components and related stress resultants, i.e., the axial force S i , shear force Q i , and bending moment M i (i = 0, 1, 2), are defined to follow the standard sign conventions designated in elementary mechanics of materials [103]. For the representative segmental element of the upper adherend layer [See Figure 2a], the static equilibrium equations in terms of stress resultants are ΣF x = 0 : ΣF y = 0 : tions in terms of stress resultants are The static equilibrium equations of the representative segmental element of the adhesive layer [See Figure 2b] are written as : The static equilibrium equations of the representative segmental element of the adhesive layer [See Figure 2b] are written as ΣF y = 0 : Figure 3a shows the configuration of an ABSLJ. A uniform shear traction t 0 is applied at the right end of the lower adherend, and the left end of the upper adherend is fixed with a distance far away from the bonding area. The geometrical and mechanical parameters are designated as: Bonding length L, thicknesses h i (i = 0, 1, 2), Young's modulus E i (i = 0, 1, 2), Poisson's ratio ν i (i = 0, 1, 2), and coefficients of thermal expansion α i (i = 0, 1, 2). Herein, subscripts 0, 1, and 2 denote the adhesive layer, upper adherend, and lower adherend, respectively. Similarly, width b of the joint is considered to be unity. The coordinate systems are adopted as follows. The x-coordinate is selected from the left end of the joint to direct along the layer axis; y 0 , y 1 , and y 2 are the vertical coordinates with the corresponding origins attached at the centroids of the cross-section of the adhesive layer, and the upper and lower adherends, respectively. It can be expected that high interfacial normal (peeling) and shear stresses can be triggered at the adherend free-ends due to the mismatch of the material properties across the adherend interfaces as illustrated in Figure 3c. Such high interfacial mechanical or thermomechanical stresses are responsible for the typical debonding failure of ABJs.  In reality, ABJs subjected to have mechanical or thermomechanical loads are typically in a general three-dimensional (3D) stress state. To simplify the modeling process, in this study the ABJs are treated in the plane-stress state and without residual stresses in the initial load-free state. Thus, the mechanical and thermomechanical stresses can be treated separately according to the method of superposition. In addition, the stress results obtained in the plane-stress state can be conveniently converted to those in the plane-strain state by simply replacing the Young's moduli

Stress Resultants of an Adhesively Bonded Single-Lap Joint (ABSLJ) Subjected to a Shear Force
It is the unique feature of the stress-function variational method to define the shear and normal (peeling) stresses at the interface between the upper adherend and the adhesive layer as two independent interfacial stress functions to be determined: Similarly, the interfacial shear and normal (peeling) stresses at the interface between the lower adherend and the adhesive layer are assumed to be another two independent interfacial stress functions to be determined: Thus, the shear-free conditions at the adherend edges at x = 0 and L stand for And Furthermore, the physical conditions of the axial and shear tractions and the bending moments at the ends of the upper and lower adherends and the adhesive layer specify the rest force conditions at the right and left ends of the three-layered joint segment corresponding to the specific ABJs. In the case of an ABSLJ subjected to a shear force at the right end of the lower adherend as shown in Figure 3, the corresponding stress resultants at the right and left ends of the three layers can be expressed as below.
In the above, t 0 and m 0 are the average shear traction and resultant bending moment per unit width in the joint. Given the interfacial shear and normal stress functions to be determined at the upper and lower surfaces of the adhesive layer (i.e., τ 1 = f 1 , τ 2 = f 2 , σ 1 = g 1 , and σ 2 = g 2 ), the axial force S 1 , shear traction Q 1 , and bending moment M 1 at an arbitrary location x of the upper adherend can expressed with BCs (13a-f) and the integration of Equations (1)-(3) as Similarly, with BCs (13g-l) and integration of Equations (4)-(6), the axial force S 2 , shear traction Q 2 , and bending moment M 2 at an arbitrary location x of the lower adherend of the ABSLJ can be determined as Furthermore, the axial force S 0 , shear traction Q 0 , and bending moment M 0 at an arbitrary location x at the adhesive layer can be determined via integrating Equations (7)-(9) with BCs (13m-r) as

Planar Stresses in the Adherends and Adhesive Layer of an ABSLJ
In the process of stress-function variational method, the procedure for determining the planar stresses in the adherends and adhesive layers is standardized [2,[23][24][25][26] such that the axial stress in each layer is assumed to be linearly varying according to the classic Euler-Bernoulli beam theory while the corresponding shear and transverse normal stresses are determined to satisfy the 2D static equilibrium equations. Thus, the axial stress of the upper adherend of the ABSLJ can be expressed as The corresponding shear stress τ (1) y 1 x in the upper adherend can be determined via integrating the 2D static equilibrium equation: with respect to y 1 from an arbitrary location y 1 to the top surface at y 1 = h 1 /2: which yields τ (1) where the traction-free BC τ (1) y 1 x (h 1 /2) = 0 is evoked. In addition, the transverse normal stress σ (1) y 1 y 1 in the upper adherend can be further determined via integrating the 2D static equilibrium equation: with respect to y 1 from an arbitrary location y 1 to the top surface at y 1 = h 1 /2 as which leads to σ (1) In the integration process above, traction BC σ (1) Similarly, the axial normal stress σ (2) xx , shear stress τ (2) y 2 x , and transverse normal stress σ (2) y 2 y 2 in the lower adherend of the ABSLJ can be determined as In derivations of (31) and (32), integrations for τ (2) y 2 x (x, y 2 ) and σ (2) y 2 y 2 (x, y 2 ) are made with the upper and lower limits as y 2 and −h 2 /2, respectively, and traction BCs of In derivations of (34) and (35), integrations for τ y 0 y 0 (x, y 0 ) are made with the upper and lower limits as y 0 and −h 0 /2, respectively, and traction BCs of (34) and (35) are automatically consistent with τ (0) in which the minus sign prior to f 1 (x) is due to the sign conversion of stress components in the theory of elasticity.
The above derivations of the stress components indicate that with the assumption of axial normal stresses varying linearly across the adherend and adhesive layers of the ABSLJ, the corresponding statically compatible shear and transverse normal stresses have piecewise parabolic and cubic distributions across these layers, respectively. In addition, such stress fields satisfy all the traction BCs at the ends of the adherend and adhesive layers and the stress continuity across the interfaces between the adherend and adhesive layers. Such a process can also be conveniently extended for determining the stress components in multi-layered ABJs.

Governing Equations of the Interfacial Stress Functions and Their Solution of an ABSLJ
With the planar stress components in the adherend and adhesive layers of the ABSLJ, the total strain energy of the ABSLJ (0 ≤ x ≤ L) can be expressed as [2] In the above, ε (i) xx and ε (i) yy (i = 0, 1, 2) are the axial and transverse normal strains of the adhesive layer and the upper and lower adherends, respectively, which are defined according to the generalized Hooke's law of an isotropic, linearly thermoelastic solid (in the plane-stress state): where α i (i = 0, 1, 2) are coefficients of thermal expansion of the adhesive layer, upper, and lower adherends, respectively, and ∆T is the uniform temperature change of the joint from the reference temperature of free thermomechanical stress state. In the present case of an ABSLJ, strain energy (36) is an energy functional with respect to the four unknown interfacial stress-functions f i (i = 1, 2) and g i (i = 1, 2) adopted above. According to theorem of minimum complimentary strain energy of an elastic body, the total strain energy (or complimentary strain energy) reaches a stationary point at static equilibrium of the present linearly elastic ABSLJ (with given tractions at free ends), which corresponds to the necessary condition in terms of variation of the strain energy (36) with respect to the four unknown stress functions [2,3,[23][24][25][26] i.e., xx δε (2) xx + δσ (2) xx ε (2) xx + σ (2) yy δε (2) yy + δσ (2) yy ε (2) yy ] + 2(1+υ 2 ) where δ is the mathematical variational operator with respect to either f i (i = 1, 2) or g i (i = 1, 2). By adopting the same notations and procedure as used in our previous studies [2,3,[23][24][25][26], plugging the stress components (23), (26), (29), and (30)- (35) and normal strains (37) and (38) into (40), and evoking the variational operation and several steps of algebraic simplification, it yields the four interfacial stress functions f 1 , f 2 , g 1 , and g 2 that need to satisfy a system of four coupled fourth-order ODEs of constant coefficients: where {Φ} 4×1 is a dimensionless interfacial stress function vector that is defined as In Equation (41) and where h 02 = h 0 /h 2 ,h 12 = h 1 /h 2 ,e 20 = E 2 /E 0 ,e 21 = E 2 /E 1 .
In addition, {D} 4×1 is a dimensionless 4×1 mechanical or thermomechanical load vector: in which It needs to be mentioned that except for expression (47), the above expressions of matrices [A], [B], and [C] carry the same expressions as determined in our previous studies [2,3], which demonstrates the high potential to universally utilize the present stressfunction variational method for high-accuracy interfacial stress analysis of ABJs. Hereafter, by following the standard procedure as formulated in [2,3], the solution to the system of governing ODEs (41) To solve the system of homogenous ODEs (48), assume that the general formal solution {Ψ} can be expressed as where λ and {Ψ 0 } are, respectively, the eigenvalue and eigenvector of a generalized eigenvalue problem corresponding to (49) such that This generalized eigenvalue problem can be further converted into a standard eigenvalue problem by introducing As a result, the generalized eigenvalue problem (53) is reduced to a standard eigenvalue problem as where I is a 4 × 4 identity. This standard eigenvalue problem can be solved numerically by using popular numerical algorithms available in the literature, e.g., the eig( ) function available in Matlab ® , etc. Consequently, the formal solution (48) can be expressed as where Ψ k 0 (k = 1, 2, . . . , 8) are eigenvectors (the first four elements of each column vector) corresponding to eigenvalues λ k (k = 1, 2, . . . , 8), respectively, and c k and d k (k = 1, 2, . . . , 8) are 16 real-valued or complex coefficients to be determined to satisfy 16 traction BCs at the ends of the adherend layers (12a), (12b), and (13a)-(13l): In addition, traction BCs (13m)-(13r) are satisfied automatically once the interfacial stress functions satisfy the set of governing equations and traction BCs (56a)-(56p). Furthermore, subjected to a shear traction of uniform force density t 0 as shown in Figure 3, differentiation of the particular solution (50) with respect to the dimensionless spatial coordinate is: Thus, substitution of (55) and (57) into the reduced BCs (56a)-(56p) leads to a set of 16 simultaneous algebraic equations for determining the unknown c k and d k (k = 1,2, . . . ,8): 0 (L/h 2 )/dξ, (58h) where Ψ k,1 0 , Ψ k,2 0 , Ψ k,3 0 , and Ψ k,4 0 (k = 1,2, . . . ,8) are respectively the first to fourth elements of the k-th eigenvector, and Φ 0 , and Φ (4) 0 are respectively the first to fourth elements of the particular solution vector {Φ 0 } as given in (50). Once the unknown coefficients c k and d k are determined, plugging (55) into (43a)-(43d) yields the shear and normal (peeling) stresses on the upper and bottom interfaces of the ABSLJ as

Model Validation by Finite Element Method (FEM)
Let us first validate the present semi-analytic stress-function variational method by using a commercially available FEM software package (ANSYS ® ) to determine and compare the interfacial normal and shear stresses of an ABSLJ subjected to a uniformly distributed shear force in the plane-stress state as shown in Figure 3a [29]. Previously, efforts of FEM-based validation were conducted by the authors in the cases of a bonded joint and an adhesively bonded single-strip joint (ABSSJ) [2,23] and by others in the literature [36,37], which indicated that this method is a high-efficiency and accurate technique for interfacial stress analysis of bonded joints and ABSSJs. In the present case, the linearly elastic ABSLJ is assumed to be made of an upper steel adherend (E 1 = 210 GPa, υ 1 = 0.293) and a lower aluminum adherend (E 2 = 70 GPa, υ 2 = 0.345), which are adhesively bonded together through an epoxy-type adhesive layer (E 0 =10 GPa, υ 0 = 0.40). The adherends and adhesive layer have the same width, and the rest geometries of the joint are: h 1 = 2.0 mm (steel), h 2 = 2.0 mm (aluminum), h 0 = 0.2 mm (adhesive), and L = 20 mm [See Figure 3a]. The uniformly distributed shear force with the force density of t 0 = 1.0 MPa is applied onto the right edge of the lower adherend. The linear elastic analysis is performed by using ANSYS ® to determine the stress field of the present ABSSJ, and no nonlinear geometrical effects such as large deformations and defections are considered, which are corresponding to the assumptions of small, linearly elastic deformations as adopted in the present stressfunction variational method. During the finite element (FE) modeling, the three-layered joint segment as shown in Figure 3d is considered. The external forces applied to the FE model are generated as follows. A pair of uniformly distributed shear-forces t 0 of unity density (t 0 =1.0 MPa) are applied on both the left edge of the upper adherend and the right free-edge of the lower adherend, and a linearly varying force is applied to the left edge of the upper adherend following the flexural stress variation of σ xx = −m 0 y 1 /E 1 , which has the bending moment −m 0 as the stress resultant. To constrain the rigid motion of the joint segment in the numerical process, both the horizontal and vertical displacements of one corner node are fixed. Besides, four-node elements (PLANE182) and mapped uniform quadrilateral meshes are utilized in the FE modeling. To comparatively study the varying trend of the singular interfacial stresses near the free edges of the ABSLJ, two refined mesh sizes (i.e., the quadrilateral elements with the dimensions of 0.1 mm × 0.1 mm and 0.05 mm × 0.05 mm, respectively) are adopted at the free-edges.
As a result, variations of the interfacial shear and normal stresses at the upper and lower surfaces of the adhesive layer with the distance from the left free-edge are plotted in Figure 4a,b. It can be observed that the shear stresses predicted by the present model exactly satisfy the shear-free condition at two free edges of the joint as enforced by this method while the shear stresses predicted by FEA do not satisfy such conditions due to the limitation of FEM. In addition, although the present method is significantly distinct from FEM, both the interfacial shear and normal stresses predicted by the present method are very close to those predicted by FEA, especially for the free-edge interfacial normal stresses. This demonstrates the high capability of the present semi-analytic method for interfacial stress analysis of ABJs. Figure 4a also indicates the noticeable difference of the interfacial shear stress across the thin adhesive layer, especially at the locus where the peak shear stress appears, which breaks the shear-lagging assumption as adopted in classic models of ABJs [21,22] and many follow-ups, while Figure 4b does not show noticeable variation of the normal stress across the thin adhesive layer. Besides, as a matter of fact, sharp corners are formed at the free edges of two adhesively bonded layers in the joint, where stress singularities exist in the view of theoretical analysis [59][60][61][62]. Thus, numerical results of the interfacial shear and normal stresses predicted by FEA are expected to be asymptotically infinite near the free edges with the meshes to be refined with decreasing dimensions. As no singularity functions are adopted in the present method, no singular stresses appear, although very high stress variations are predicted.
Furthermore, due to existence of the bending moment and shear force at the left edge of the upper adherend of the ABSLJ as shown in Figure 3a, the peak values of both the interfacial shear and normal stresses are detected near the left free edge. Such a high interfacial normal (peeling) stress is responsible mainly for the potential debonding failure of the ABSLJ; meanwhile, the high interfacial shear stress further enhances the effective failure stress, i.e., von Mises stress. Therefore, the present stress-function variational method can also be used to validate other analytic and numerical methods for interfacial stress analysis of ABSLJs, and this method can be further used for scaling analysis for design and failure analysis of ABSLJs. of the ABSLJ; meanwhile, the high interfacial shear stress further enhances the effec failure stress, i.e., von Mises stress. Therefore, the present stress-function variatio method can also be used to validate other analytic and numerical methods for interfa stress analysis of ABSLJs, and this method can be further used for scaling analysis design and failure analysis of ABSLJs.

Scaling Analysis of Interfacial Shear and Normal Stresses of ABSLJs due to Mechanical Loads
With the high-efficiency, accurate semi-analytic stress-function variational met developed above, it is convenient to examine the effects of elastic properties and geo tries of the adherend and adhesive layers on the interfacial shear and normal stresse the ABSLJ. Below the effects of thickness and Young's modulus of the adhesive laye the interfacial shear and normal stresses in the present ABSLJ are investigated. For c venience of the scaling analysis, the adherend length and Young's modulus ratios fixed as L/h2 = 8 and E2/E1 = 1/3 (approximately equal to the ratio of aluminum to ste and Poisson's ratios of the upper and lower adherends are fixed as υ1 = 0.293 (steel) υ2 = 0.345 (aluminum), respectively. For the adhesive layer, four thickness ratios (h0/ 0.1, 0.25, 0.5, and 1.0) and two Young's modulus ratios (E0/E2 = 1/10 and 1/5) are adop and the Poisson's ratio is fixed at υ0 = 0.4 (thermosetting epoxy) in the entire scaling a ysis. The joint is treated in the plane-stress state. Figures 5 and 6 show variations of normalized interfacial shear stress τ/t0 and normal (peeling) stress σ/t0 at the upper lower surfaces of the adhesive layer with the dimensionless distance x/h2 from the le

Scaling Analysis of Interfacial Shear and Normal Stresses of ABSLJs Due to Mechanical Loads
With the high-efficiency, accurate semi-analytic stress-function variational method developed above, it is convenient to examine the effects of elastic properties and geometries of the adherend and adhesive layers on the interfacial shear and normal stresses of the ABSLJ. Below the effects of thickness and Young's modulus of the adhesive layer on the interfacial shear and normal stresses in the present ABSLJ are investigated. For convenience of the scaling analysis, the adherend length and Young's modulus ratios are fixed as L/h 2 = 8 and E 2 /E 1 = 1/3 (approximately equal to the ratio of aluminum to steel), and Poisson's ratios of the upper and lower adherends are fixed as υ 1 = 0.293 (steel) and υ 2 = 0.345 (aluminum), respectively. For the adhesive layer, four thickness ratios (h 0 /h 2 = 0.1, 0.25, 0.5, and 1.0) and two Young's modulus ratios (E 0 /E 2 = 1/10 and 1/5) are adopted, and the Poisson's ratio is fixed at υ 0 = 0.4 (thermosetting epoxy) in the entire scaling analysis. The joint is treated in the plane-stress state. Figures 5 and 6 show variations of the normalized interfacial shear stress τ/t 0 and normal (peeling) stress σ/t 0 at the upper and lower surfaces of the adhesive layer with the dimensionless distance x/h 2 from the left to the right adherend end at a fixed adherend thickness ratio (h 1 /h 2 = 1.0), four adhesive thickness ratios (h 0 /h 2 = 0.1, 0.25, 0.5 and 1.0), and two adhesive Young's modulus ratios (E 0 /E 2 = 1/10 and 1/5), respectively. A related Matlab ® code for the results shown in Figure 6 is attached in Appendix A.
From Figures 5 and 6, the numerical results of all the cases show that the shear stresses satisfy the shear-free condition at each adherend end. Due to the unique geometries of the ABSLJ as shown in Figure 3a, the peak values of both the interfacial shear and normal stresses appear at the left adherend end of the joint. In addition, at a fixed modulus ratio E 0 /E 2 , the peak values of the interfacial shear and normal stresses increase with decreasing thickness ratio of the adhesive layer, i.e., h 0 /h 2 . For instance, when the thickness ratio h 0 /h 2 decreases from 1.0 to 0.1, the peak value of the interfacial shear stress is doubled, i.e., a thinner adhesive layer corresponds to the higher interfacial shear stresses. Such an observation is obvious since the adhesive layer is much more compliant than the two adherends in this study and in broad ABJs used in engineering practices. Therefore, the adhesive layer is able to offer sufficient large deformations to suppress the larger deformations in the adherends. Furthermore, within the range of the adhesive-layer thickness ratios from 0.1 to 1.0 in this study, the peak value of the interfacial normal (peeling) stress at the left adherend end behaves with a more complicated feature as follows.
On one side, the peak value of the interfacial normal stress at the left free-end increases significantly with decreasing h 0 /h 2 , i.e., an obvious size effect of the adhesive layer exists on the interfacial normal stress variation. On the other side, the peak value of the interfacial normal stress has a significant variation across the adhesive layer, which again indicates that the classic shear-lagging assumption of the adhesive layer and ignorance of the normal stress variation across the adhesive layer are oversimplified in establishing accurate mechanics models of ABJs for stress and strength analysis. Moreover, Figures 5 and 6 also demonstrate the noticeable effect of the adhesive layer modulus on both the interfacial shear and normal stresses. The higher the adhesive layer modules, the higher are the interfacial shear and normal stresses. Thus, detailed scaling analysis of the interfacial stresses of ABSLJs can be applicable for rational design, structural optimization, and failure analysis of related ABJs. With the present semi-analytic method, such scaling analysis becomes more convenient by comparison with purely computational stress analysis of ABJs based on FEA. Consequently, it needs also to be mentioned that choice of the above adhesive thickness ratio h 0 /h 2 for scaling analysis can be in a broader range as in many cases the thickness of adhesive layers can be much thinner compared to that of the adherend layers. In the numerical process, when h 0 /h 2 is relatively small, say < 1%, the condition number of the coefficient matrix of the set of 16 final linear algebraic equations can be large, which may lead to some additional numerical errors. As a numerical example, it was verified that the present semi-analytic method and related numerical algorithm as shown in the attached Matlab ® code in Appendix A are capable of predicting the accurate interfacial stresses for ABJs with the adhesive thickness ratio h 0 /h 2 as small as 0.1%.

Interfacial Shear and Normal Stress Analysis of an Adhesively Single-Sided Joint (ASSJ) Subjected to Uniaxial Tension
The above effective semi-analytic stress-function variational method can be further extended for determining the stress fields, in particular the interfacial shear and normal stresses, in other broad types of ABJs subjected to mechanical and thermomechanical loads. Below the interfacial shear and normal (peeling) stress analysis of an ASSJ subjected to uniaxial tension and uniform temperature change are considered, respectively. The solving procedure in this case is very similar to that described in Sections 2.1, 2.2 and 2.4. For the purpose of simplification, herein just a list of the main results for an ASSJ is shown in Figure 7a. The ASSJ under consideration is assumed to be made of two rectangular elastic adherend layers adhesively bonded through an elastic adhesive layer. All the symbols, coordinates, and material properties are defined following those adopted for the above ABSLJ in Sections 2.1 and 2.2. Thus, by following the problem-solving procedure for an ABSLJ, the stress resultants, planar stress components, and related traction BCs of the ASSJ subjected to uniaxial tension in the lower adherent layer can be expressed in the following (in plane-stress state). ABSLJ in Sections 2.1-2.2. Thus, by following the problem-solving procedure for an A SLJ, the stress resultants, planar stress components, and related traction BCs of the AS subjected to uniaxial tension in the lower adherent layer can be expressed in the followi (in plane-stress state).
The axial force S2(x), shear force Q2(x), and bending moment M2(x) in the lower a herend of the ASSJ are The axial force S 2 (x), shear force Q 2 (x), and bending moment M 2 (x) in the lower adherend of the ASSJ are The axial force S 0 (x), shear force Q 0 (x), and bending moment M 0 (x) in the adhesive layer of the ASSJ are σ (1) The planar axial normal stress σ (2) xx , shear stress τ (2) xy 2 , and lateral normal stress σ (2) y 2 y 2 in the lower adherend of the ASSJ are The planar axial normal stress σ Again, by evoking the principle of minimum complimentary strain energy of the ASSJ as demonstrated in Section 2.4, it again results in the governing ODE (41) with the same coefficient matrices A, B, and C as given in (44), while the matrix D in (47) is modified as the following. {D} 4×1 is a dimensionless 4 × 1 mechanical or thermomechanical load vector: with The particular solution to the governing ODE (43) and its derivative in this case of ASSJ are Correspondingly, the traction BCs for the governing ODE (41) of the current ASSJ can be expressed as The rest procedure for solving the interfacial stress functions for an ASSJ is the same as that for an ABSLJ as aforementioned. Thus, the numerical solutions for the interfacial shear and normal stresses in an ASSJ can be conveniently determined via modifying the relevant terms in the Matlab ® code as given in Appendix A, which again demonstrates the high efficiency and universality of the stress-function variational method for interfacial stress analysis of various types of ABJs.

Model Validation and Numerical Examples of Interfacial Shear and Normal Stresses of an ASSJ Subjected to Uniaxial Tension Model Validation by Finite Element Method (FEM)
Similar to Section 2.5.1, the FEM software package ANSYS ® is further utilized for validating the stress-function variational method for ASSJs. Detailed computational processes based on ANSYS ® for the present ASSJ are the same as demonstrated in Section 2.5.1. The upper and lower adherends are assumed to be made from steel and aluminum, respectively, and the adhesive layer is assumed to be made from thermosetting epoxy. The elastic properties of each layer and the mesh size as well as the element type adopted in FEA are the same as those used in Section 2.5.1. In addition, in the present case of an ASSJ under uniaxial tension as shown in Figure 7a, both the joint geometries and mechanical loads are symmetrical. Thus, only the right half symmetric joint is needed for validation by FEA, and single-point displacement constraint is applied to constrain the lateral rigid motion of the joint. Figure 8 shows the detailed numerical results of the interfacial shear and normal (peeling) stresses at the upper and lower surfaces of the adhesive layer of the ASSJ that are determined by using the present semi-analytic as well as FEA with two mesh sizes. Similar to those of ABSLJ in Section 2.5.1, FEA-based validation endorses high accuracy and high efficiency of the present stress-function variational method for interfacial stress analysis of ABJs.

Scaling Analysis of Thermomechanical Interfacial Shear and Normal Stresses in an AS Subjected to Uniform Temperature Change
Thermomechanical stress analysis of an ASSJ [See Figure 7a] due purely to a unifo temperature change is equivalent to that of an adhesively bonded thermostat. Thermom chanical stress analysis of bonded and adhesively bonded thermostats (e.g., aluminu molybdenum) has been extensively studied by a number of investigators [2,10,11,1 Hereafter, a simple scaling analysis is conducted for examining the effect of the adhes thickness on the thermomechanical interfacial stress variations in an ASSJ. During scaling analysis based on the present method, the geometrical configuration and ela properties of the above ASSJ are adopted as follows: Material properties: Steel: E1 = 2 GPa, υ1 = 0.293, α1 = 11.3 × 10 −6 /°C (coefficient of thermal expansion), h1 = 2.0 mm; Alu num: E2 = 70 GPa, υ2 = 0.345, α2 = 23.6 × 10 −6 /°C, h2 = 2.0 mm; Epoxy-type adhesive: E0 = GPa, υ0 = 0.4, α0 = 57.6 × 10 −6 /°C; length of bonding line: L = 50 mm. The three layers of ASSJ are assumed to carry the same width and in the plane-stress state, and it is furt assumed that there are no residual stresses at the reference temperature. Figure 9 sho the interfacial shear and normal stresses at the upper and lower surfaces of the adhes layer of the ASSJ subjected purely to a uniform temperature increase of ΔT = 120 o from reference temperature at four thickness ratios of the adhesive layer to the lower adhere layer, i.e., h0/h2 = 0.25, 0.5, 1.0, and 2.0. From Figure 9, it can be clearly observed that Scaling Analysis of Thermomechanical Interfacial Shear and Normal Stresses in an ASSJ Subjected to Uniform Temperature Change Thermomechanical stress analysis of an ASSJ [See Figure 7a] due purely to a uniform temperature change is equivalent to that of an adhesively bonded thermostat. Thermomechanical stress analysis of bonded and adhesively bonded thermostats (e.g., aluminummolybdenum) has been extensively studied by a number of investigators [2,10,11,14]. Hereafter, a simple scaling analysis is conducted for examining the effect of the adhesive thickness on the thermomechanical interfacial stress variations in an ASSJ. During the scaling analysis based on the present method, the geometrical configuration and elastic properties of the above ASSJ are adopted as follows: Material properties: Steel: E 1 = 210 GPa, υ 1 = 0.293, α 1 = 11.3 × 10 −6 / • C (coefficient of thermal expansion), h 1 = 2.0 mm; Aluminum: E 2 = 70 GPa, υ 2 = 0.345, α 2 = 23.6 × 10 −6 / • C, h 2 = 2.0 mm; Epoxy-type adhesive: E 0 = 10 GPa, υ 0 = 0.4, α 0 = 57.6 × 10 −6 / • C; length of bonding line: L = 50 mm. The three layers of the ASSJ are assumed to carry the same width and in the plane-stress state, and it is further assumed that there are no residual stresses at the reference temperature. Figure 9 shows the interfacial shear and normal stresses at the upper and lower surfaces of the adhesive layer of the ASSJ subjected purely to a uniform temperature increase of ∆T = 120 • from the reference temperature at four thickness ratios of the adhesive layer to the lower adherend layer, i.e., h 0 /h 2 = 0.25, 0.5, 1.0, and 2.0. From Figure 9, it can be clearly observed that the thermomechanical shear stresses satisfy the shear-free conditions at the two free edges of the ASSJ. Both high interfacial shear and normal (peeling) stresses are highly localized near the free-edges, where debonding failure most likely happens, and these interfacial shear and normal stresses decay rapidly at the locations off the free edges. In addition, in each case, both the interfacial shear and normal stresses increase with increasing thickness ratio h 0 /h 2 of the adhesive layer, i.e., a relatively thicker adhesive layer corresponds to higher stress concentrations near the free edges. This is due to the fact that the adhesive layer carries a much higher coefficient of thermal expansion. At a fixed temperature change ∆T, a thicker adhesive layer yields higher mismatch deformations between the adhesive layer and either the upper or lower adherend layer of the ASSJ, i.e., higher mismatch strains and high interfacial stresses. Furthermore, the peak value of the peeling stress is nearly double that of the interfacial shear stress at the free edges, thus this thermomechanical peeling stress is responsible mainly for the debonding failure. In microelectronics packaging, the high mismatch of the coefficients of thermal expansion between chips and binding polymers commonly leads to the high thermomechanical stresses, which are responsible for the thermomechanical fatigue failure of microelectronic systems. It needs to also be mentioned that as the present analysis is based on a linearly elastic ASSJ model, all predicted interfacial shear and normal stress only depend on the thickness and modulus ratios. Therefore, the present stress-function variational method also demonstrates the great potential for efficient and accurate thermomechanical interfacial stress analysis of ABJs and adhesively bonded thermostats. ASSJ model, all predicted interfacial shear and normal stress only depend on the thickness and modulus ratios. Therefore, the present stress-function variational method also demonstrates the great potential for efficient and accurate thermomechanical interfacial stress analysis of ABJs and adhesively bonded thermostats.

Conclusions
A generalized stress-function variational method was successfully formulated for stress analysis of ABJs subjected to mechanical and thermomechanical loads. The key of the present method for determining the interfacial shear and normal stresses in an ABJ

Conclusions
A generalized stress-function variational method was successfully formulated for stress analysis of ABJs subjected to mechanical and thermomechanical loads. The key of the present method for determining the interfacial shear and normal stresses in an ABJ was to introduce two unknown shear and normal stress functions at each interface of the ABJ, and then all the stress components in the adherends and the adhesive layer of the ABJ could be expressed in terms of these unknown interfacial stress functions within the framework of the classic Euler-Bernoulli beam theory and linear elasticity. By evoking the principle of minimum complimentary strain energy of the ABJ, a set of governing ODEs were obtained, which are well-conditioned with proper traction BCs and can be solved conveniently by using eigenfunction method via designing a compact, universal Matlab ® code. Two examples of ABSLJs and ASSJ were demonstrated to show high efficiency and accuracy of the present semi-analytic method. It can be concluded that the set of governing ODEs formulated by the authors in this and previous studies is universal for all the three-layered ABJs as shown in Figure 1, except for matrix D given in (47), which corresponds to the specific force BCs and specific type of ABJ. Consequently, though the present work only considers two simple three-layered ABJs, the present semianalytic method can be conveniently extended for analysis of layer-by-layer interfacial stresses in multi-layered ABJs and ABJs made of composite laminates and smart materials, fragmentation of fibers and surface coatings, and mechanical durability of layered materials and flexible electronics, etc.