Acessibilidade / Reportar erro

Elastic and Viscoelastic Analysis of Stress Distribution in Single Lap Joints under Tensile Loading: Analytical and Numerical Approaches

Abstract

In this study, stress variations in cylindrical single lap joints under tensile loading were investigated as both elastic and viscoelastic(time-dependent) solutions. Analytical modeling utilized the Pugno-Carpinteri method and the Lubkin-Reissner method, while numerical modeling employed the Abaqus analysis program. The viscoelastic solution was derived from the associated elastic solution application using the Alfrey correspondence principle. Transformed elastic modulus values for the viscoelastic solution were obtained through methodologies proposed by Schapery and the generalized Maxwell-Wiechert model. The analyses revealed that maximum stress occurred at the initiation of loading, with stress values decreasing over time at the adhesive edges, and an increase in stiffness observed in the central region of the adhesive. In single lap joints, it was found that overall in-plane stiffness remained constant independent of the adhesive modulus value.

Keywords
Cylindrical bonding joint; Stress analysis; Finite element analysis; Least-squares method; Alfrey’s correspondence principle

Graphical Abstract

1 INTRODUCTION

Adhesive bonding joints represent a primary method for fastening structural elements, driven by the growing demand for lightweight, high-quality, and cost-effective products. This joining system ensures load-bearing capability while preserving structural integrity. Adhesives find extensive applications across industries such as aviation, automotive, shipping, and power plants. They offer comparable or superior joint strength compared to traditional bonding methods like screws, rivets, and welding. Despite the ability of adhesives to distribute loads and stresses uniformly over the bond surface, thereby mitigating stress concentrations at specific points, a degree of stress concentration within the adhesive layer persists (Aimmanee, 2018Aimmanee, S. (2018). A unified analysis of adhesive bonded cylindrical coupler joints. Applied adhesive bonding in science and technology, 11-31.). Cylindrical adhesive joints constitute a focal point of investigation within the literature.

Viscoelasticity entails the study of materials whose properties change over time. The advancement of viscoelastic analysis has been further propelled by the emergence of “synthetic polymer materials” exhibiting specific characteristics. While polymer materials behave elastically at low temperatures and high deformation rates, they exhibit viscous behavior at high temperatures and low deformation rates. Additionally, investigations into viscoelasticity hold significance not only in mechanics but also in biomechanics, as many biomaterials, such as heart tissue, muscle tissue, and cartilage, display a viscoelastic response. In stress analysis of cylindrical bonded adhesives, the adhesive is typically considered as an elastic material. However, in reality, most polymers used as adhesives are viscoelastic materials with stress-strain properties that change over time. Alfrey's correspondence principle provides a method for utilizing solutions known in the context of elasticity to address viscoelastic problems.

In the literature, some examples of elastic and viscoelastic studies related to stress analysis of cylindrically bonded adhesives are presented below.

In the field of analitic analysis studies, Lubkin and Reissner conducted groundbreaking research on cylindrical joints. They provided a closed-form solution for peel and shear stresses in pipe joints subjected to axial loading. Opting to utilize prescribed boundary conditions instead of relying on initial values within the standard Laplace transformation framework, they thus established their innovative methodology as a cornerstone in this area. Additionally, Adams and Peppiart (1977)Adams, R., & Peppiart, N. (1977). Stress analysis of adhesive bonded tubular lap joints. J. Adhesion, Vol.9, 1-18. along with Dragoni and Goglio (2013)Dragoni, E., & Goglio, L. (2013). Adhesive stresses in axially loaded tubular bonded joints-part-I. International journal of adhesion and adhesives, Vol.47, 35-45. delved into the intricacies of shear and peel stresses within the adhesive through detailed finite element analyses. Emphasizing the predominance of shear and peel stresses over other stress components, they indicated the negligible nature of the latter. Concurrently, by affirming the efficacy of Lubkin and Reissner's model in accurately predicting both peel and shear stresses, they validated the significance of their approach. Furthermore, Gunay et al. (1998)Gunay, D., Aydemir, A., & Özer, H. (1998). Stress analysis in cylindrical adhesive lap joint of pressure vessels. 8th international machine design and production conference semptember 9-11, Ankara, 339-349. conducted numerical examinations on stresses in cylindrical adhesive joints, enriching the understanding of local stress phenomena by extensively investigating stress concentration levels at adhesive edges using the finite element method. Meanwhile, Saraç (2021)Saraç, İ. (2021). Boru yapıştırma bağlantılarında farklı tasarım parametrelerinin yapıştırıcı tabakasında gerilme dağılımına etkisinin sayısal olarak araştırılması. Uludağ Üniversitesi Mühendislik Fakültesi Dergisi, 26(1). conducted a comprehensive numerical exploration by examining various design parameters (particularly overlap length, pipe wall thickness, and adhesive thickness) that elucidate the complex interaction between structural geometry and stress distribution in adhesive joints, providing valuable insights into stress distribution within the adhesive region. Expanding on this theme, Campilho et al. (2022)Campilho, R., Pinheiro, A., Moreira, R., & Sanchez-Arce, I. (2022). Validation of theoretical models for the strengh prediction of tubular adhesive joints. Procedia Structure Integrity, 41:60-71. performed a comparative analysis using three different adhesives with varying properties Araldite AV138, Araldite 2015, and SikaForce 7752 on cylindrical adhesive joints, highlighting the occurrence of highest stresses in the overlap region due to the brittle nature of Araldite AV138. Khan and Kumar (2017)Khan, M., & Kumar, S. (2017). Mitigating edge effects in adhesively bonded composite tubular lap joints under axial, pressure and torsional loading via stiffness-tailoring. 21st International conference on composite materials, 20-25th August. Xi'an. emphasized the relatively superior load transfer behavior of adhesive pipe joints compared to their traditional counterparts, with maximum stress concentrations typically occurring at the ends of the overlap region, underscoring pragmatic considerations in adhesive joint design and analysis. Additionally, Kumar and Pandey observed that the “static strength” of joints with hybrid adhesive joint lines was higher than that of joints with single lap joint lines. With parametric studies, higher bond strength was also achieved for the optimum two-adhesive bond line ratio. It has also been analytically demonstrated that the overall in-plane stiffness of the joint remains unaffected by the modulus rating of the bonding line adhesive (Kumar & Pandey, 2010Kumar, S., & Pandey, P. (2010). Behaviour of bi-adhesive joints. (24).).

In the field of viscoelastic studies, Alwar and Nagaraja were the first to consider the time-dependent properties of the adhesive by assuming its viscoelastic nature. Concerning stress analysis, they compared finite element method results obtained using Prony series to represent the relaxation modulus for two different adhesives with the Lubkin-Reissner theory. Additionally, in other articles, they compared finite element method findings for flat-lap joints bonded with adhesive with the Goland-Reissner theory (Nagaraja & Alwar, 1976Nagaraja, Y., & Alwar, R. (1976). Viscoelastic analysis of an adhesive tubular joint. Journal of Adhesion, Vol.8, 76-92.; Alwar & Nagaraja, 1980Alwar, R., & Nagaraja, Y. (1980). Viscoelastic analysis of an adhesive bonded lap joints. Computers & Structures, Vol.11, 621-627.). Delale and Edoğan (1981)Delale, F., & Edoğan, F. (1981). Viscoelastic analysis of adhesively bonded joints. Journal of Applied Mechanical, Vol. 48, 331-338. conducted analyses assuming that the pipes were elastic and the adhesive was linearly viscoelastic. After formulating the problem, they utilized the standard Laplace transformation method for solving a special scenario that combines two identical pipes with a three-parameter viscoelastic solid adhesive. They calculated stress distribution in the adhesive layer under three different external loads: membrane loading, bending, and transverse shear loading. Their results indicated that normal stress (peel) in the adhesive was lower and degraded slower than shear stress. Furthermore, they mentioned that in the case of a viscoelastic adhesive, the elastic response could also be determined for t = +0 and t = ∞ situations using limit theorems for inverse Laplace transformations in the s domain. Yilmazyurt et al. (2022)Yilmazyurt M., Eyüpreisoğlu, S., Okyar, A. F., & Namlı, O. C. (2022). Viscoelastic characterization and mechanical hysteresis of commercial grade polypropylene. Journal of Polytechnic. utilized stress relaxation tests to understand how the viscoelastic behavior of polypropylene material could be mathematically modeled. They compared the time-dependent stress-strain response of polypropylene material through finite element solution in Marc software with the analytical solution of mathematical models. Cheng et al. (2013)Cheng, F., Özsoy, Ö., & Reddy, N. (2013). Finite element modeling of viscoelastic behaviour and interface damage in adhesively bonded joints. Scriver Publishing LLC, 23-46. performed failure analysis of bonded joints under tensile load for both elastic and viscoelastic adhesives at specific loading rates. They mentioned that compared to an elastic adhesive, a viscoelastic adhesive could increase the stiffness of the adhesive-bonded system. Additionally, based on the cohesive zone model results for adhesive-pipe interfaces, they observed that higher loading rates accelerated crack initiation and propagation at adhesive-pipe interfaces. Abouel-Kasem (2014)Abouel-Kasem, A. (2014). Analysis and design of viscoelastic adhesively bonded tubular joint. Journal of Enginering and Apllied Sciences, Vol.1, 13-23. conducted stress analysis under quasi-static internal pressure considering the viscoelastic properties for assemblies of pipes bonded with six different adhesives. They examined the effects of connection geometry and loading conditions on fatigue and equivalent stresses. He noted that initially applied high stress values could cause edge cracks in the adhesive layer. Shishesaz and Reza (2013)Shishesaz, M., & Reza, A. (2013). The effect of viscoelasticity of polymeric adhesives on shear stress distribution in a single lap joint. The Journal of Adhesion, 89: 859-880. investigated the effect of adhesive viscoelasticity on the shear stress distribution in the adhesive layer of a single-lap joint under tensile load. They found a relationship between analytical results for shear stress in the adhesive layer and finite element results and mentioned that the ratio of adhesive's viscous effect to shear modulus had an inverse effect on the peak shear stress developed in the joint. Tapar and Reza (2021)studiedTapar, M. H., & Reza, A. (2021). Long- term shear stress distribution in adhesively bonded tubular joints under tensile load using the time- temperature superposition principle. The Journal of Adhesion, Vol. 97, 328-345. the distribution of long-term shear stress due to temperature in assemblies of pipes bonded with adhesive under tensile load. They obtained the temperature-dependent relaxation modulus of the viscoelastic adhesive using the time-temperature superposition principle. They noted that maximum shear stress occurred at low temperatures. Özer and Çay (2022)Özer, H., & Çay, M. (2022). Viscoelastic model to evaluate the shear stresses in the bi-adhesively bonded lap joints. The Brazilian Society of Mechanical Sciences and Engineering, Vol. 44, 518. developed a viscoelastic model to modify the viscoelastic behavior of hybrid adhesives and examine their effects on time-dependent stress distributions. They stated that the viscoelastic property in hybrid adhesive assemblies was more effective in flexible adhesives at the ends and less effective in stiff adhesive in the middle.

In this study, stress distributions along the cylindrical adhesive overlap joint subjected to tensile load were obtained analytically and numerically using Alfrey's correspondence principle. In elastic solutions, the problem was analytically modeled employing the Pugno-Carpinteri method (Pugno & Carpinteri, 2003Pugno, N., & Carpinteri, A. (2003). Tubular adhesive joints under axial load. Journal of Applied Mechanics, ASME, Vol. 70, 832-839.) and the Lubkin and Reissner method (Lubkin & Reissner, 1956). Numerical analyses were conducted using Abaqus CAE. In viscoelastic analyses, the Prony series model, also known as the Wiechert model, was employed. Elastic and time-dependent variations of stresses along the overlap length of cylindrical bonding joints were examined analytically and numerically, and the results were compared.

2 ANALYSİS MODEL AND ALFREY’S CORRESPONDENCE PRINCIPLE

2.1 Analysis model: Geometry, Materials and Boundary Conditions

Two different adhesive situations were examined using flexible and rigid adhesives in the tubular single lap joint shown in Figure 1. The inner radius of the inner tube (ri=8 mm), the thickness of both tubes and the adhesive (h1,2=2 mm, ha=0.2 mm), and the length of each adhesive (L=50 mm) were kept constant throughout the analyses. For tubular adhesive joints with homogeneous adhesive, the geometric specification was taken from (Campilho, Pinheiro, Moreira, & Sanchez-Arce, 2022) and the material properties were taken from (Özer & Çay, 2022). Applied load: 7615.208 N.

Figure 1
Cylindrical lap joint model and boundary conditions.
Figure 2
Elementary free body diagrams for the joint.

Where F denotes the axial force, E1, E2 and Ea represent the elastic moduli, G1, G2, and Ga denote the shear moduli, r1, r2, and a are the radii of the middle surface, υ1, υ2, and υa are the Poisson’s ratios of the adherends and adhesive, and L=2c is the overlap length. The right end of the outer tube is constrained in all directions (Ux=Ur=UƟ=0), while the left end of the inner tube is constrained only in the radial direction (Ur=0). The cylindrical coordinate system is denoted by (x, r, θ).

2.2 Alfrey’s Correspondence Principle

The concept of exploiting elastic solutions to find viscoelastic solutions was first proposed by Alfrey Turner in 1944 and later generalized in 1950 and 1955. This principle was further developed by W. T. Read and E. H. Lee respectively. In viscoelastic stress analysis, an equivalent problem in elastic theory is obtained by eliminating the time variable through the Laplace transform parameter. Then, by substituting a real-time elastic modulus into the solution of this elastic problem, the viscoelastic solution of the problem is obtained (Turner, 1948Turner, A. J. (1948). Mechanical Behaviour of High Polymers. New York: Intescience Publishers.; Lee, 1954Lee, E. (1954). Stress analysisi in viskoelastic bodies. Brown University, 183-190.).

In the Alfrey correspondence principle, it is stipulated that the viscoelastic equations and boundary conditions, when transformed into the Laplace domain with zero initial condition, must maintain formal equivalence with the equations governing an elastic system possessing identical geometry. Additionally, while applied loads and displacements may vary with time, the bonding surface or interface must remain unaltered throughout. Consequently, the Laplace transform operates solely on the time variable (t) and remains constant with respect to the spatial variable (x). Broadly speaking, this methodology can be applied to any scenario wherein load functions, including boundary conditions, can be expressed as a product of spatial and temporal components. It is asserted that employing this approach facilitates the utilization of the extensive literature available in elasticity theory for tackling viscoelastic problems.

Alfrey’s correspondence principle comprises four key steps:

In the initial step, the elastic solution of the problem is determined. Subsequently, Laplace transforms are substituted for elastic variables in the stress equations. The transformed modulus value is obtained by multiplying the Prony series Laplace transform by “s”.

Finally, upon transitioning from the Laplace domain to the time domain, the obtained solution represents the viscoelastic solution of the problem.

3 ANALYSIS METHODS

3.1 Analytical Analysis methods

3.1.1 Pugno- Carpinteri Analytical Method:

The Pugno-Carpinteri method involves the analytical calculation of shear stresses within the bonding region, specifically along the midpoint of the adhesive thickness, in co-axial pipe bonding joints (see Figure 1). This method aims to obtain stress distributions by primarily considering pure shear while neglecting bending stresses.

In the analytical model, the equation of shear stress is expressed by the following equations. Taken together with other formulas provided by the reference (Pugno & Carpinteri, 2003)

τrxx=-F2πRC1αeαx-C2αe-αx,(1)

Where, F: Tensile load; ti: Thickness; R, a, ri:radii. The constants C1 and C2 are obtained under boundary conditions and are defined together with other constants as follows.

C1=e-αce-2αc - e2αc+βeαc - e-αce-2αc - e2αc ,(2)
C2=eαce2αc - e-2αc+βe-αc - eαce2αc - e-2αc ,(3)
α=KE1A1+E2A2E1A1E2A2 ,(4)
β=E1A1E1A1+E2A2 ,(5)
K=2π RGata ,(6)

From classical strength to average shear stress

τm=F2π*R*L ,(7)

In this study, Stress Concentration Factor (SCF) analysis was also conducted. SCF is defined below:

SCF=τrxτm=- F2πR αC1eαx-C2e-αxF4π R c=-2αcC1eαx-C2e-αx,(8)

3.1.2 Lubkin- Reissner Analytical Method:

In the Lubkin and Reissner method considered the bending moment in tubular lap joints under axial tensile load (see Figure 1). It provided solutions for both shear stress (τrx) and normal stresses (σr) throughout the thickness of the adhesive layer. Outer tube (adhered 2) and inner tube (adhered 2).

Where, Axial (F1, F2), transverse (Q1, Q2), hoop (N1, N2) forces and bending moments (M1, M2), per unit length in Fig. 3 and Fig. 4. Adhesive shear (τrx) and peeling (σrr) stress balance equations were written. The following formulas are taken from reference (Lubkin & Reissner, 1956).

Figure 3
Basic free body diagrams for cylindrical joint, (r: Radial, θ: Tangential, x: Axial).
Figure 4
Basic free body diagrams for cylindrical joint, (r: Radial, θ: Tangential, x: Axial).

Differential Equation Set

Using the unknown terms v1, v2 and F0 defined below, three sets of simultaneous differential equations are obtained.

v1=F2πa1-υ12E12ct12g1z , v2=F2πa1-υ22E22ct22g2z, F0=F2πag3z (9)

v1, v2 and F0 given in Equation 15 are unknown (corresponding to the normalized vertical/radial displacements and the difference of axial forces in the pipes, respectively).

Three sets of simultaneous differential equations, two fourth-order and one second-order, and 10 boundary conditions are obtained by mathematical operations of three dimensionless functions g1(z), g2(z), and g3(z); Here z=((x+c))/2c is a dimensionless abscissa.

Then, shear and peeling stress distribution graphs were found by substituting three-dimensional g1(z), g2(z), and g3(z) into the equations. The obtained stress distributions represent the stress states in the middle thickness of the adhesive.

g 1 i v + K 1 4 + γ 11 4 g 1 - γ 12 4 g 2 - 3 a a 1 g 3 ' ' - 3 Λ 1 a a 1 g 3 = - 3 Λ 1 a a 1 (10)
g 2 i v + K 2 4 + γ 22 4 g 2 - γ 21 4 g 1 - 3 a a 2 g 3 ' ' + 3 Λ 2 a a 2 g 3 = - 3 Λ 2 a a 2 (11)
g 3 ' ' - a B 1 2 a 1 + a B 2 2 a 2 g 3 - B 2 2 g 2 ' ' + B 1 2 g 1 ' ' + Λ 2 B 2 2 g 2 - Λ 1 B 1 2 g 1 = a B 2 2 a 2 - a B 1 2 a 2 (12)
Λ i = 2 υ i 2 c / t i 2 t i / a i i = 1 , 2 Λ i = 2 υ i 2 c / t i 2 t i / a i i = 1 , 2 K i 4 = 12 ( 1 - υ i 2 ) 2 c / t i 4 t i / a i 2 i = 1 , 2 B i 2 = ( 1 - υ i 2 ) 2 c / t i 2 t i G a / t a E i i = 1 , 2 γ i j 4 = 12 ( 1 - υ j 2 ) 2 c / t j 3 2 c / t i a E a t j / a i E j t a i , j = 1 , 2 (13)

Boundary conditions for the set of differential equations

At the left end of the overlap, (z=0), (Lubkin & Reissner, 1956Lubkin, J., & Reissner, E. (1956). Stress distribution and design data for adhesive lap joints between circular tubes. ASME, 1213-1221.):

F2=Q2=M2=0, v1=V12λ13D1-M12λ12D1-υ1F2πE1t1 and dv1dx=Q12λ12D1-M1λ1D1
g 3 0 = - 1
g 2 0 ' ' = 0
g 2 0 ' ' ' - 3 a a 2 g 3 0 ' = 0 (14)
g 1 0 ' ' ' - 2 K 1 g 1 0 ' ' + K 1 2 g 1 0 ' - 3 a a 1 g 3 0 ' = 0
g 1 0 ' ' - 2 K 1 g 1 0 ' + K 1 2 g 1 0 = - 2 3 υ 1 a a 1 1 - υ 1 2

At the right end of the overlap, (z=1), (Lubkin & Reissner, 1956Lubkin, J., & Reissner, E. (1956). Stress distribution and design data for adhesive lap joints between circular tubes. ASME, 1213-1221.):

F1=Q1=M1=0, v2=-Q22λ23D2-M22λ22D2-υ2F2πE2t2 ve dv2dx=Q22λ22D2+M2λ2D2
g 3 1 = 1
g 1 1 ' ' = 0
g 1 1 ' ' ' - 3 a a 1 g 3 1 ' = 0 (15)
g 2 1 ' ' ' + 2 K 2 g 2 1 ' ' + K 2 2 g 2 1 ' - 3 a a 2 g 3 1 ' = 0
g 2 1 ' ' + 2 K 2 g 2 1 ' + K 2 2 g 2 1 = - 2 3 υ 2 a a 2 1 - υ 2 2

Laplace Transform and Solution of Differential Equations

By substituting the geometric and material values given in Figure 1, the differential equations are derived and then transformed using the Laplace transform. These equations are subsequently expressed in matrix form.

s 4 G 1 s - s 3 g 1 0 - s 2 g 1 ' 0 - s g 1 ' ' 0 - g 1 ' ' ' 0 + K 1 4 + γ 11 4 G 1 s - γ 12 4 G 2 s - 3 a a 1 s 2 G 3 s - s g 3 0 - g 3 ' ' 0 - 3 Λ 1 a a 1 G 3 s = - 1 s 3 Λ 1 a a 1
s 4 + K 1 4 + γ 11 4 G 1 s - γ 12 4 G 2 s - 3 a a 1 s 2 + Λ 1 G 3 s = 1 s s 4 g 1 0 + s 3 g 1 ' 0 + s 2 ( g 1 ' ' 0 - 3 a a 1 g 3 0 ) + s ( g 1 ' ' ' 0 - 3 a a 1 g 3 ' 0 ) - 3 Λ 1 a a 1
s 4 G 2 s - s 3 g 2 0 - s 2 g 2 ' 0 - s g 2 ' ' 0 - g 2 ' ' ' 0 + K 2 4 + γ 22 4 G 2 s - γ 21 4 G 1 s - 3 a a 2 s 2 G 3 s - s g 3 0 - g 3 ' ' 0 - 3 Λ 2 a a 2 G 3 s = - 1 s 3 Λ 2 a a 2
- γ 21 4 G 1 s + s 4 + K 2 4 + γ 22 4 G 2 s - 3 a a 2 s 2 - Λ 2 G 3 s = 1 s s 4 g 2 0 + s 3 g 2 ' 0 + s 2 ( g 2 ' ' 0 - 3 a a 2 g 3 0 ) + s ( g 2 ' ' ' 0 - 3 a a 2 g 3 ' 0 ) - 3 Λ 1 a a 1
s 2 G 3 s - s g 3 0 - g 3 ' 0 - B 2 2 a 2 + B 1 2 a 1 a G 3 s - ( B 2 2 ( s 2 G 2 s - s g 2 0 - g 2 ' 0 + B 1 2 s 2 G 1 s - s g 1 0 - g 1 ' 0 ) + Λ 2 B 2 2 G 2 s - Λ 1 B 1 2 G 1 s = 1 s B 2 2 a / a 2 - B 1 2 a / a 1
- B 1 2 s 2 + Λ 1 G 1 s - B 2 2 s 2 - Λ 2 G 2 s + s 2 - B 2 2 a 2 + B 1 2 a 1 a G 3 s = 1 s s 2 g 3 0 - B 1 2 g 1 0 - B 2 2 g 2 0 + s g 3 ' 0 - B 1 2 g 1 ' 0 - B 2 2 g 2 ' 0 - B 1 2 a a 1 - B 2 2 a a 2
A 3 x 3 * G ( s ) 3 x 1 = B 3 x 1
s 4 + K 1 4 + γ 11 4 - γ 12 4 - 3 a a 1 s 2 + Λ 1 - γ 21 4 s 4 + K 2 4 + γ 22 4 - 3 a a 1 s 2 + Λ 1 - B 1 2 s 2 + Λ 1 - B 2 2 s 2 - Λ 2 s 2 - B 1 2 a 1 + B 2 2 a 2 a G 1 s G 2 s G 3 s = 1 s s 4 g 1 0 + s 3 g 1 0 ' + s 2 g 1 0 ' ' - 3 a a 1 g 3 0 + s g 1 0 ' ' ' - 3 a a 1 g 3 0 ' - 3 Λ 1 a a 1 s 4 g 2 0 + s 3 g 2 0 ' + s 2 g 2 0 ' ' - 3 a a 2 g 3 0 + s g 2 0 ' ' ' - 3 a a 2 g 3 0 ' - 3 Λ 2 a a 2 s 2 g 3 0 - B 1 2 g 1 0 - B 2 2 g 2 0 + s g 3 0 ' - B 1 2 g 1 0 ' - B 2 2 g 2 0 ' - B 1 2 a a 1 - B 2 2 a a 2

The equation of the tenth degree is subject to five boundary conditions specified at z=0. The remaining five boundary conditions are determined using the five equations provided at z=1. (See Appendix A Appendix A clc; clear all; % Analytical Shear and Peel Stress Distributions under tensile load % Geometric and materials properties format short E1=207000; E2=207000;Ea=2700;ta=0.2;L=20;F=1269.203; R1i=8; R1o=10;R2i=10.2; R2o=12.2; v1=0.3;v2=0.3;va=0.367; t1=2;t2=2; Ga=Ea/(2*(1+va)); % Lubkin Reisner Method syms g_1_0 g_1_0_t g_2_0 g_2_0_t g_3_0_t syms t s % Inner pipe properties (Steel). L=20; c=L/2; a1=(R1i+R1o)/2 ;t1=2; E1=207000; v1=0.3; % Outer pipe properties(Steel). L=20; c=L/2; a2=(R2i+R2o)/2 ;t2=2; E2=207000; v2=0.3; % Adhesive properties. a=10.2;ta=0.2; Ea=2700; va=0.367; F=1269.203; tau=1; Ga=Ea/(2*(1+va)); K1=(12*(1-v1^2)*(2*c/t1)^4*(t1/a1)^2)^(1/4); K2=(12*(1-v2^2)*(2*c/t2)^4*(t2/a2)^2)^(1/4); B1=((1-v1^2)*(2*c/t1)^2*((t1*Ga)/(ta*E1)))^(1/2); B2=((1-v2^2)*(2*c/t2)^2*((t2*Ga)/(ta*E2)))^(1/2); g11=(12*(1-v1^2)*(2*c/t1)^3*(2*c/t1)*((a*Ea*t1)/(a1*E1*ta)))^(1/4); % g11=gamma11; g22=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t2)*((a*Ea*t2)/(a2*E2*ta)))^(1/4); % g22=gamma22; g12=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t1)*((a*Ea*t2)/(a1*E2*ta)))^(1/4); % g12=gamma12; g21=(12*(1-v1^2)*(2*c/t1) ^3*(2*c/t2)*((a*Ea*t1)/(a2*E1*ta)))^(1/4); % g21=gamma21; L1=2*v1*(2*c/t1)^2*(t1/a1); % L1=lambda1 L2=2*v2*(2*c/t2)^2*(t2/a2); % L2=lambda2 g30=-1; g20_2t=0; % fixed values at z=0 % T1=g10^”-3*a/a1*g30; T1=sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2))-3*(a/a1)*(g30); % T2=g10^”’-3*(a/a1)*g30^’ T2=sqrt(2)*K1*(sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2)))-(K1^2)*g_1_0_t; % T3=g20^”-3*(a/a2)*g30; T3=g20_2t-3*(a/a2)*g30; % T4=g20^”’-3*(a/a2)* g30^’=0; T4=0; % T5=g30-(B1^2)*g_1_0-(B2^2)*g_2_0; T5= g30-(B1^2)*g_1_0-(B2^2)*g_2_0; % T6=g30^’-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t; T6=g_3_0_t-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t; A=[s^4+K1^4+g11^4,-g12^4,-3*(a/a1)*(s^2+L1);-g21^4,s^4+K2^4+g22^4,-3*(a/a2)*(s^2-L2);-(B1^2)*(s^2+L1), -(B2^2)*(s^2-L2),s^2-(a/a1)*(B1^2)-(a/a2)*B2^2]; B=(1/s)*[g_1_0*s^4+g_1_0_t*s^3+T1*s^2+T2*s-3*L1*(a/a1);... g_2_0*s^4+g_2_0_t*s^3+T3*s^2+T4*s-3*L2*(a/a2);... T5*s^2+T6*s-(a/a1)*(B1^2)+(a/a2)*B2^2]; G=A\B; g=ilaplace(G); disp(g) g3z=g(3); g1z=g(1); g2z=g(2); % Two function and three function derivatives unknown % By arranging the g functions according to the 5 boundary conditions given in z=1, a 5x5 matris equation is formed. [M]*{D}=[P] M=zeros(5); q=subs(g3z,t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(1,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(1,1)=1-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g1z,t,2),t,1); q=subs(diff(g1z,t,2),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(2,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(2,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g1z,t,3),t,1)-3*(a/a1)*subs(diff(g3z,t),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(3,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(3,1)=-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g2z,t,3),t,1)+sqrt(2)*K2*subs(diff(g2z,t,2),t,1)+(K2^2)*subs(diff(g2z,t),t,1)-3*(a/a2)*subs(diff(g3z,t),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(4,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(4,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g2z,t,2),t,1)+sqrt(2)*K2*subs(diff(g2z,t),t,1)+(K2^2)*subs(g2z,t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(5,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(5,1)=(-2*sqrt(3)*v2*((a/a2)*(1/sqrt(1-v2^2))))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); % g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t values % The D values are substituted in the g functions and the actual g1z, g2z, g3z functions are found. D=M\P; z=linspace(0,1,20); g=subs(g,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},D'); % D' column vector g1z=g(1); g2z=g(2); g3z=g(3); figure(1) % Plotting g1. for i=1:length(z) g1(i)=double(subs(g1z,t,z(i))); end plot(z, g1,'-k','linewidth',2) title(' g_1(z) Dimensionless function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_1(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_1(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','best'); text(0.1,-0.019,'(a) ','Fontsize',10, 'Fontweight', 'bold'); grid minor figure(2) % Plotting g2. for i=1:length(z) g2(i)=double(subs(g2z,t,z(i))); end plot(z, g2,'-r','linewidth',2) title(' g_2(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_2(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_2(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast'); grid minor text(0.1,-0.018,'(b) ','Fontsize',10, 'Fontweight', 'bold'); figure(3) % Plotting g3. for i=1:length(z) g3(i)=double(subs(g3z,t,z(i))); end plot(z, g3,'-k','linewidth',2) title(' g_3(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_3(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_3(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast'); grid minor text(0.1,-0.8,'(c) ','Fontsize',10, 'Fontweight', 'bold'); figure(4) % Plotting Shear Stress. for i=1:length(z) g3(i)=double(subs(g3z,t,z(i))); shear(i)=double(subs(diff(g3z,t),t,z(i)))/2; end plot(z, shear,'-r','linewidth',2) title('\tau_{rx}- Shear Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b'); xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('\tau_{rx}/ \tau_m Normalized Stress ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k') legend('\tau_{rx} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold'); xlim([0 1]) text(0.3,1.45,[('Lubkin/Reissner = '), num2str(max(shear),4)],'Fontsize',10, 'Fontweight', 'bold') text(0.1,0.7,'(d) ','Fontsize',12, 'Fontweight', 'bold'); grid minor figure(5) % Plotting Peel Stress. for i=1:length(z) g2(i)=double(subs(g2z,t,z(i))); v_2(i)=tau*((1-v2^2)/(E2*t2^2))*(2*c)^3*g2(i); g1(i)=double(subs(g1z,t,z(i))); v_1(i)=tau*((1-v1^2)/(E1*t1^2))*(2*c)^3*g1(i); peel(i)=(Ea/ta)*(v_2(i)- v_1(i)); end plot(z, peel,'-r','linewidth',2) title('\sigma_{rr} - Peel Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b'); xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('\sigma_{rr} / \tau_m Normalized Peel Stress', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k') legend('\sigma_{rr} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold'); text(0.3,0.75,[('Lubkin/Reissner = '), num2str(max(peel),4)],'Fontsize',10, 'Fontweight', 'bold') text(0.1,-0.3,'(e) ','Fontsize',12, 'Fontweight', 'bold'); grid minor ).

[M]*{D}=[P] D=M-1P=g1og1o'g2og2o'g3o'5X1= -0.0062 0.0164 0.0031-0.0858 3.51665X1

By substituting the D values into B, the vector Gi(s)3x1​ is first obtained. Subsequently, the inverse Laplace transform is applied to derive the dimensionless functions gi(z)3x1.

Shear stress (τrx) and Peel stress (σrr):

Shear stress and Peel stress equations are completed with adhesive stress-strain relationships:

τ r x = G a γ a = G a t a u 2 , i n - u 1 , o u (16)

The dimensionless functions, derived from material and geometric properties, determine that the ratio of the displacement difference (∆v) to the adhesive thickness provides the radial strain. Multiplying the strain by the elastic modulus of the adhesive yields the radial stress, known as the peel stress.

σ r r = E a ε a = E a v t a = E a t a v 2 - v 1 (17)

These equations are derived from Lamé theory for thick-walled cylinders. When the thickness of the cylinders is small, the Lubkin-Reissner model becomes applicable, providing a more accurate description of the stress and strain distribution in such scenarios. When the wall thickness of the cylinders is small, the Lubkin-Reissner model becomes more appropriate and valid, as it accurately analyzes the behavior of thin-walled structures under small deformations. For thick-walled cylinders, however, the Lamé theory provides more accurate results. Therefore, when the wall thickness of the cylinders is small, more precise analyses can be achieved using the Lubkin-Reissner model.

Normalized Shear and Peel Stresses:

τ m = F 4 π * a * c = 1 (18)
τ r x τ m = F / 8 π a c g 3 ' ( z ) F / 4 π a c = g 3 ' ( z ) 2 (19)
σ r r τ m = E a G a B 2 2 2 c t 2 g 2 - B 1 2 2 c t 1 g 1 (20)
B22=1-υ222ct22t2GataE2, B12=1-υ122ct12t1GataE1(21)

The solution to the set of Differential Equations is based on the Laplace transform. The shear and peel stress distributions along the overlap length, as depicted in Figure 1, were analyzed using the Lubkin Reissner method implemented in Matlab software (Mathwork Inc. 2023Mathwork Inc. (2023). The student edition off Matlab.), employing Equations (9) through (21). The material of the pipes was uniformly defined as S235JR, and Sika Power 4720 was utilized as the adhesive. The geometric, material, and mechanical properties are presented in Table 1 and Table 2, respectively. An applied tensile load of 1269.203 N was considered for the analysis. (See Appendix A Appendix A clc; clear all; % Analytical Shear and Peel Stress Distributions under tensile load % Geometric and materials properties format short E1=207000; E2=207000;Ea=2700;ta=0.2;L=20;F=1269.203; R1i=8; R1o=10;R2i=10.2; R2o=12.2; v1=0.3;v2=0.3;va=0.367; t1=2;t2=2; Ga=Ea/(2*(1+va)); % Lubkin Reisner Method syms g_1_0 g_1_0_t g_2_0 g_2_0_t g_3_0_t syms t s % Inner pipe properties (Steel). L=20; c=L/2; a1=(R1i+R1o)/2 ;t1=2; E1=207000; v1=0.3; % Outer pipe properties(Steel). L=20; c=L/2; a2=(R2i+R2o)/2 ;t2=2; E2=207000; v2=0.3; % Adhesive properties. a=10.2;ta=0.2; Ea=2700; va=0.367; F=1269.203; tau=1; Ga=Ea/(2*(1+va)); K1=(12*(1-v1^2)*(2*c/t1)^4*(t1/a1)^2)^(1/4); K2=(12*(1-v2^2)*(2*c/t2)^4*(t2/a2)^2)^(1/4); B1=((1-v1^2)*(2*c/t1)^2*((t1*Ga)/(ta*E1)))^(1/2); B2=((1-v2^2)*(2*c/t2)^2*((t2*Ga)/(ta*E2)))^(1/2); g11=(12*(1-v1^2)*(2*c/t1)^3*(2*c/t1)*((a*Ea*t1)/(a1*E1*ta)))^(1/4); % g11=gamma11; g22=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t2)*((a*Ea*t2)/(a2*E2*ta)))^(1/4); % g22=gamma22; g12=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t1)*((a*Ea*t2)/(a1*E2*ta)))^(1/4); % g12=gamma12; g21=(12*(1-v1^2)*(2*c/t1) ^3*(2*c/t2)*((a*Ea*t1)/(a2*E1*ta)))^(1/4); % g21=gamma21; L1=2*v1*(2*c/t1)^2*(t1/a1); % L1=lambda1 L2=2*v2*(2*c/t2)^2*(t2/a2); % L2=lambda2 g30=-1; g20_2t=0; % fixed values at z=0 % T1=g10^”-3*a/a1*g30; T1=sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2))-3*(a/a1)*(g30); % T2=g10^”’-3*(a/a1)*g30^’ T2=sqrt(2)*K1*(sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2)))-(K1^2)*g_1_0_t; % T3=g20^”-3*(a/a2)*g30; T3=g20_2t-3*(a/a2)*g30; % T4=g20^”’-3*(a/a2)* g30^’=0; T4=0; % T5=g30-(B1^2)*g_1_0-(B2^2)*g_2_0; T5= g30-(B1^2)*g_1_0-(B2^2)*g_2_0; % T6=g30^’-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t; T6=g_3_0_t-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t; A=[s^4+K1^4+g11^4,-g12^4,-3*(a/a1)*(s^2+L1);-g21^4,s^4+K2^4+g22^4,-3*(a/a2)*(s^2-L2);-(B1^2)*(s^2+L1), -(B2^2)*(s^2-L2),s^2-(a/a1)*(B1^2)-(a/a2)*B2^2]; B=(1/s)*[g_1_0*s^4+g_1_0_t*s^3+T1*s^2+T2*s-3*L1*(a/a1);... g_2_0*s^4+g_2_0_t*s^3+T3*s^2+T4*s-3*L2*(a/a2);... T5*s^2+T6*s-(a/a1)*(B1^2)+(a/a2)*B2^2]; G=A\B; g=ilaplace(G); disp(g) g3z=g(3); g1z=g(1); g2z=g(2); % Two function and three function derivatives unknown % By arranging the g functions according to the 5 boundary conditions given in z=1, a 5x5 matris equation is formed. [M]*{D}=[P] M=zeros(5); q=subs(g3z,t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(1,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(1,1)=1-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g1z,t,2),t,1); q=subs(diff(g1z,t,2),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(2,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(2,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g1z,t,3),t,1)-3*(a/a1)*subs(diff(g3z,t),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(3,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(3,1)=-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g2z,t,3),t,1)+sqrt(2)*K2*subs(diff(g2z,t,2),t,1)+(K2^2)*subs(diff(g2z,t),t,1)-3*(a/a2)*subs(diff(g3z,t),t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(4,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(4,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); q=subs(diff(g2z,t,2),t,1)+sqrt(2)*K2*subs(diff(g2z,t),t,1)+(K2^2)*subs(g2z,t,1); for i=1:5 W=zeros(1,5); W(i)=1; M(5,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); end P(5,1)=(-2*sqrt(3)*v2*((a/a2)*(1/sqrt(1-v2^2))))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5))); % g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t values % The D values are substituted in the g functions and the actual g1z, g2z, g3z functions are found. D=M\P; z=linspace(0,1,20); g=subs(g,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},D'); % D' column vector g1z=g(1); g2z=g(2); g3z=g(3); figure(1) % Plotting g1. for i=1:length(z) g1(i)=double(subs(g1z,t,z(i))); end plot(z, g1,'-k','linewidth',2) title(' g_1(z) Dimensionless function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_1(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_1(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','best'); text(0.1,-0.019,'(a) ','Fontsize',10, 'Fontweight', 'bold'); grid minor figure(2) % Plotting g2. for i=1:length(z) g2(i)=double(subs(g2z,t,z(i))); end plot(z, g2,'-r','linewidth',2) title(' g_2(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_2(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_2(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast'); grid minor text(0.1,-0.018,'(b) ','Fontsize',10, 'Fontweight', 'bold'); figure(3) % Plotting g3. for i=1:length(z) g3(i)=double(subs(g3z,t,z(i))); end plot(z, g3,'-k','linewidth',2) title(' g_3(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b') xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('g_3(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); legend('g_3(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast'); grid minor text(0.1,-0.8,'(c) ','Fontsize',10, 'Fontweight', 'bold'); figure(4) % Plotting Shear Stress. for i=1:length(z) g3(i)=double(subs(g3z,t,z(i))); shear(i)=double(subs(diff(g3z,t),t,z(i)))/2; end plot(z, shear,'-r','linewidth',2) title('\tau_{rx}- Shear Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b'); xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('\tau_{rx}/ \tau_m Normalized Stress ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k') legend('\tau_{rx} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold'); xlim([0 1]) text(0.3,1.45,[('Lubkin/Reissner = '), num2str(max(shear),4)],'Fontsize',10, 'Fontweight', 'bold') text(0.1,0.7,'(d) ','Fontsize',12, 'Fontweight', 'bold'); grid minor figure(5) % Plotting Peel Stress. for i=1:length(z) g2(i)=double(subs(g2z,t,z(i))); v_2(i)=tau*((1-v2^2)/(E2*t2^2))*(2*c)^3*g2(i); g1(i)=double(subs(g1z,t,z(i))); v_1(i)=tau*((1-v1^2)/(E1*t1^2))*(2*c)^3*g1(i); peel(i)=(Ea/ta)*(v_2(i)- v_1(i)); end plot(z, peel,'-r','linewidth',2) title('\sigma_{rr} - Peel Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b'); xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); ylabel('\sigma_{rr} / \tau_m Normalized Peel Stress', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k') legend('\sigma_{rr} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold'); text(0.3,0.75,[('Lubkin/Reissner = '), num2str(max(peel),4)],'Fontsize',10, 'Fontweight', 'bold') text(0.1,-0.3,'(e) ','Fontsize',12, 'Fontweight', 'bold'); grid minor ).

Table 1
Geometric properties (Campilho, Pinheiro, Moreira, & Sanchez-Arce, 2022Campilho, R., Pinheiro, A., Moreira, R., & Sanchez-Arce, I. (2022). Validation of theoretical models for the strengh prediction of tubular adhesive joints. Procedia Structure Integrity, 41:60-71.).
Table 2
Material mechanical properties (Özer & Çay, 2022).

As a result, the elastic solution of the problem, which represents the first step of the Alfrey correspondence principle, was determined using two different approaches.

3.2 Viscoelastic Analysis

3.2.1 Generalized Maxwell Model (Wiechert Model):

In the Generalized Maxwell model depicted in Figure 5, the material comprises n discrete Maxwell elements along with a spring

Figure 5
Generalized Maxwell Model (Wiechert Model).

In the case of constant strain, the total stress is the summation of stresses within the spring and each Maxwell element. Among the available models for solids, the Wiechert model, also referred to as the Prony series model, employs as many spring-damper elements as required. Deriving the differential equation using a single Maxwell element fails to yield accurate results for real materials. A range of 5 to 15 or more Maxwell elements may be necessary to accurately capture the viscoelastic behavior of the material. Typically, “n” logarithmic time intervals are chosen across the transition range of the stress relaxation test curve, and a τk value is determined within each of these “n” intervals. Subsequently, by plotting these coordinate points at “n” positions along the experimental curve, the elastic modulus function can be readily determined using the Least Squares Method (Nagaraja & Alwar, 1976Nagaraja, Y., & Alwar, R. (1976). Viscoelastic analysis of an adhesive tubular joint. Journal of Adhesion, Vol.8, 76-92.; Brinson & Brinson, 2008Brinson, H., & Brinson, C. (2008). Polymer Engineering Science and Viscoelasticity An Introduction. New York: Springer Science-Business Media, LLC.). Thus, it is imperative to utilize a sufficient number of spring-dashpot components to precisely represent the relaxation modulus data acquired from stress relaxation testing (Machiraju, Phan, Pearsal, & Madanagopal, 2006Machiraju, C., Phan, V., Pearsal, A., & Madanagopal, S. (2006). Viscoelastic studies of human subscapularis tendon: Relaxation test and a Wiechert model. Computer Methods and Programs in Biomedicine. Computer Methods and Programs in Biomedicine, 29–33.; Yılmazyurt M., Eyüpreisoğlu, Okyar, & Namlı, 2022).

σ t = σ 1 + σ 2 + σ n + σ (22)
σ t = ε 0 E t σ t = ε 0 i = 1 n E i e - t τ i + E
σt=0=ε0E1+E2++E, σt==ε0E
Et=σtε0=E+i=1nEi e- tτi , ε0=1 için, σt=Et
E t = E + i = 1 n E i . e - 1 τ i t
E = E 0 - i = 1 n E i E t = ( E 0 - i = 1 n E i ) + i = 1 n E i . e - 1 τ i t
Et=E0-(i=1nEi-i=1nEi.e-1τi t), Et=E0-i=1nEi(1-e-1τi t)
Et=E+i=1nEi e- tτi=E0 [1-i=1ngi (1-e-1τi t)] (Prony Series)(23)
gi=EiE0, τi=ηiEi , Prony series coefficients (gi, τi)(24)

The Laplace transform of the Prony series can be easily obtained from Equation 29.

L E t = L E + i = 1 n E i e - t τ i E - s = E s + E 1 s + τ 1 + E 2 s + τ 2 + E n s + τ n (25)

In the generalized Maxwell model, the order of the differential equation increases by one for each additional Maxwell element. The differential equation is then simplified to a standard algebraic problem and solved using the Laplace transform. When the strain history is defined, the generalized Maxwell model becomes more applicable and offers improved conformity with real materials. Hooke's law is employed in the general equation depicting the viscoelastic behavior of the material.

Pσ=Qε (Generalized Hooke's law)
P=p0+p1ddt++pndndtn=i=0npididti, Q=q0+q1ddt++qndndtn=i=0nqididti
i = 0 n p i d i σ d t i = i = 0 n q i d i ε d t i (26)
p n d n σ d t n + p n - 1 d n - 1 σ d t n - 1 + + p 1 d σ d t + p 0 σ = q n d n ε d t n + q n - 1 d n - 1 ε d t n - 1 + + q 1 d ε d t + q 0 ε

In Figure 8, the arm (n) lacks a damper (q0=0), representing the single-element generalized Maxwell equation (Example, n=1)

Figure 8
(a) Plane, (b) 3D Finite element model. Loading and boundary conditions.
P 1 d σ d t + P 0 σ = q 1 ε ˙ P 1 σ ˙ + P 0 σ = q 1 ε ˙
Pσ=Qε laplace transform P-σ-=Q-ε- E-=σ-ε-=sQ-P- (27)

3.2.2 Boltzman Integral and Convolution Integral:

The Boltzmann superposition principle is instrumental in representing polymer behavior through integral equations. The Boltzmann integral, alternatively referred to as the Duhamel integral or the Hereditary integral. In the Boltzmann integral within the time domain, one of the functions is expressed in terms of (τ), while the other is expressed in terms of (t - τ). Consequently, the stress variation from inception to equilibrium state of the viscoelastic element can be observed, contingent upon its time history. The convolution integral facilitates the transformation from the time domain to the Laplace domain, thereby simplifying problem-solving, and possesses the commutative property within the convolution integral.

σ t = 0 - t E t - τ d ε τ d τ d τ (28)

τ represents any instant of time, where H(t) denotes the direct delta function of the Heaviside (Step) function and its first derivative, δ(τ).

σt=0-t Et-τdετdτdτ, ετ=ε0 Htdετdτ=ε˙=ε0 δτ , δτ=1
σ t = 0 - t E t - τ d ε τ d τ d τ = 0 - t E t - τ ε 0 δ τ d τ = ε 0 0 - t E t - τ d τ
Shifting Property: Et-τδtdτ=Et
σ t = ε 0 0 - t E t - τ d τ = ε 0 E t (29)
Convolution integral ft*gt=gt*ft=L0tfτgt-τdu=f-sg-s
L f ' t * g t = 0 t d f u d u g t - u d u = s f - s g - s
σ t = 0 - t E t - τ ε ˙ d τ = E - s s ε - s = s E - s ε - s = E - * s ε - s
σ t = ε 0 E t σ - s = E - * s ε - s
E = s E - s = E - * s (30)

Substituting this “E” value for the elastic modulus in the elastic solution enabled the derivation of the viscoelastic solution pertinent to the application.

3.3. Transformed Adhesive Elastic Modulus

R. A. Schapery has indicated that when the elastic solution is known either analytically or numerically, the viscoelastic solution can be easily computed using the actual time-dependent elastic properties of the material, regardless of the complexity of the elastic solution. Schapery presented two methods for utilizing the “transformed elastic modulus” within the framework of Alfrey's correspondence principle. In order to determine the transformed elastic modulus E-*s Schapery proposed two methods (Schapery, 1961Schapery, R. (1961). Two simple approximate methods of laplace transform inversion for viscoelastic stress analysis. Pasadena, California: Aeronautical Research Laboratory Wright Air Development Division.).

It has been noted that in the case of viscoelastic adhesive, the elastic response t = +0 and t =∞ can also be determined using limit theorems for inversion of Laplace transforms in the s-domain (Arifoğlu, 2021Arifoğlu, U. (2021). Diferansiyel Denklemler ve Matlab Çözümleri. Istanbul: Alfa BasımYayım.). It has been observed that in the context of viscoelastic adhesives, the elastic response at t = 0 and t = ∞ can be ascertained by employing limit theorems for Laplace transform inversion in the s-domain (Delale & Edoğan, 1981; Arifoğlu, 2021).

3.3.1. Direct Method:

A mathematical property of the Laplace transform was utilized, with the Laplace parameter set to s=0,5. The Prony series in Equation 29 was derived from either the stress relaxation test or the viscoelastic values provided in the table. The transformed modulus E-*s was obtained by multiplying “s” by the Laplace transform of the Prony series.

E - s = E s + i = 1 n E i s + E i η i = E + i = 1 n s * τ i * E i 1 + s * τ i
E - * s = s * E - s = E + i = 1 n s * τ i * E i 1 + s * τ i E + i = 1 n 0,5 * τ i * E i t + 0,5 * τ i (31)

3.3.2 Least-Squares Method:

It is not as simple as the direct method, but it has advantages;

Firstly, the derivative is not confined to functions that change slowly with logarithmic time.

Secondly, the time dependence is represented by a simple exponential series, which can be readily applied in the Duhamel superposition integral to compute responses to prescribed loads and displacements without requiring step functions.

Thirdly, the accuracy of the inversion can be enhanced by incorporating additional terms into the series. This is evident from Equation 23.

The viscoelastic properties of two different adhesives in Table 3 were used in both analytical and digital analyses. When the Equation 23 is written for Table 3

Table 3
Viscoelastic material properties (Özer & Çay, 2022Özer, H., & Çay, M. (2022). Viscoelastic model to evaluate the shear stresses in the bi-adhesively bonded lap joints. The Brazilian Society of Mechanical Sciences and Engineering, Vol. 44, 518.), (Karlsson, 2014Karlsson, P. (2014). Determination of viscoelastic properties of adhesives. Linnaeus University, Sweden: Master’s Thesis in Structural Engineering.)
E t = E + i 3 E i e x p - t τ i
E r e l τ 0 = E + E 1 + E 2 + E 3
E r e l τ 1 = E + E 1 e - τ 1 / τ 1 + E 2 e - τ 1 / τ 2 + E 3 e - τ 1 / τ 3
E r e l τ 2 = E + E 1 e - τ 2 / τ 1 + E 2 e - τ 2 / τ 2 + E 3 e - τ 2 / τ 3 (32)
E r e l τ 3 = E + E 1 e - τ 3 / τ 1 + E 2 e - τ 3 / τ 2 + E 3 e - τ 3 / τ 3
1 1 1 1 1 e - τ 1 / τ 1 e - τ 1 / τ 2 e - τ 1 / τ 3 1 e - τ 2 / τ 1 e - τ 2 / τ 2 e - τ 2 / τ 3 1 e - τ 3 / τ 1 e - τ 3 / τ 2 e - τ 3 / τ 3 4 x 4 E E 1 E 2 E 3 4 x 1 = E r e l 0 E r e l τ 1 E r e l τ 2 E r e l τ 3 4 x 1

3.3.3. Limit theorems: Initial value and final value:

The limit theorem involves the utilization of Laplace transform, elastic function, and its first derivative.

Initial value theorem:

L E ' t = s * E s - E 0 s * E s = L E ' t + E 0
limss*Es=limsLE't+E0 ,
lim s s * E s = lim s 0 e - s t E ' t d t + E 0 = 0 lim s e - s t E ' t d t + E 0
(limse-st=0), lims00E'tdt=0 limss*Es=0+E0
lim s s * E s = E 0 lim s s * E s = lim t 0 E t = E 0 (33)

Final value theorem:

s * E s = L E ' t + E 0
lim s 0 s * E s = lim s 0 L E ' t + E 0
lims0s*Es=lims00e-stE'tdt +E0=0lims0e-stE'tdt+E0,
(lims0e-st=1), lims0s*Es=0E'tdt+E0=Et0+E0
lim s 0 s * E s = E - E 0 + E 0 = E
lim s 0 s * E s = lim t E t = E (34)

3.4. Numerical Modeling:

Numerical investigations were carried out using the finite element method. Numerical analyzes were performed using Abaqus software and plane and three-dimensional (3D) modeling.

Results obtained from finite element analyzes were compared with Pugno-Carpinteri and Lubkin- Reissner analytical approaches.

Plane and Three-Dimensional Modeling:

In the plane (axisymmetric) finite element model, a four-node axisymmetric element was employed. The finite element model consists of 30.000 elements. A detailed view of the mesh near the overlap tips in the plane finite element mesh structure is provided in Figure 7a. Since it is modeled axisymmetrically, the boundary conditions are represented in the polar coordinate system.

Figure 7
(a) Plane, (b) 3D. Mesh structure of the joint in finite element models.

In the 3D finite element model, the finite element mesh was created using hexahedral elements. The 3D model comprises 42.000 elements and 50.837 nodes, utilizing reduced integration. The 3D finite element mesh structure is shown in Figure 7b. The results of the 3D finite element analysis are converted from Cartesian coordinates to cylindrical coordinates.

Load and Boundary Conditions:

One end of the outer tube is constrained in all directions (Ux=Ur=UƟ=0), while the inner tube is fixed in the radial direction from the inside of its end. Axial tensile load was applied in plane modeling and tensile pressure was applied in 3D modeling, and the loading and boundary conditions are seen in Figure 8a, b.

4 RESULTS AND DISCUSSIONS

The shear stress distributions along the overlap length shown in Figure 1 were obtained by applying two analytical and two numerical methodologies. The obtained shear stress values were normalized by scaling according to the average shear stress value. While the material forming the pipes is uniform and defined as S235JR (steel), Sika Power 4720 was used as the adhesive.

The geometric characteristics of the pipes, as well as the mechanical attributes of both the pipes and the adhesive employed in the modeling procedure, are delineated in Table 1 and Table 2, correspondingly. Figures 9a and 9b portray the shear stress distributions obtained from the Pugno-Carpinteri and Lubkin-Reissner methodologies, alongside 2D and 3D finite element analyses. The analytical approach proves inadequate in capturing the stress diminution at the adhesive terminations. A convergence between analytical and numerical planar outcomes becomes evident as one approaches the adhesive terminations. As illustrated in Figure 9, a discernible incongruity emerges between the findings of the planar and 3D numerical simulations at the adhesive terminations. This incongruity is believed to arise from the three-dimensional stress state present at these specific junctures. They are not necessarily more accurate as they are affected by various other numerical factors involved in modeling the material interface.

Figure 9
a, b The Pugno-Carpinteri and Lubkin-Reissner methods, in addition to 2D and 3D finite element analysis.

The application of Equation 28-36, which shows instantaneous and long-term stress changes under the condition (ε=1) for Sika Power 4720 adhesive, was carried out based on the data presented in Figure 6 and Table 3. The differential equation of the model and its derivatives at t=0 were determined. Then, stress relaxation graphs were created for the differential equation and derivative values, Prony series and limit theorem.

E - s = E s + E 1 s + E 1 η 1 + E 2 s + E 2 η 2 + E 3 s + E 3 η 3 E - s = s Q - P - = σ - ε - = s 2700 s 3 + 3500 s 2 + 250.2 s + 0.5597 s 4 + 1.446 s 3 + 0.1163 s 2 + 0.0003 s
s 4 σ - + 1.446 s 3 σ - + 0.1163 s 2 σ - + 0.0003 s σ - = s 2700 s 3 ε - + 3500 s 2 ε - + 250.2 s ε - + 0.5597
(P3s3σ-=P3d3σdt3=P3σ'''), (Terms can be written in three ways).
σ 4 ' + 1.446 σ ' ' ' + 0.1163 σ ' ' + 0.000292 σ ' = 2700 ε 4 + 3500 ε ' ' ' + 250.2 ε ' ' + 0.5597 ε '

Since ε = Constant, all derivatives are zero (In Stress Relaxation)

Figure 6
Generalized Maxwell Model (Seven-element).

In Figure 10a, calculations were conducted using Matlab software, and three graphs were generated. All three graphs exhibit identical characteristics. Figure 10b illustrates that the relaxation modulus attains equilibrium over time.

Figure 10
a, b Application for Sika Power 4720.
E 4 ' + 1.446 E ' ' ' + 0.1163 E ' ' + 0.000292 E ' = 0

E (0) = 2700

E ' 0 = - 364.29825826571072
E ' ' 0 = 463.0544372640083
E ' ' ' 0 = - 627.4251866571625

In Figures 11a and 11b, the function graphs corresponding to the relaxation modules at four different time points for Sika Power 4720 and Sika Fast 5215 adhesives were generated using Equations 31 and 32. The coefficients of the elastic function are closely aligned with the values presented in Table 3.

Figure 11
a, b The relaxation function corresponding to 4 coordinate points derived from the stress relaxation test.
E = E - * s = s * E - s E + 0,5 * τ 1 * E 1 t + 0,5 * τ 1 + 0,5 * τ 2 * E 2 t + 0,5 * τ 2 + 0,5 * τ 3 * E 3 t + 0,5 * τ 3
E P o w e r t = 1916.3417 + 246.7985 * e - t τ 1 + 301.2980 * e - t τ 2 + 235.5616 * e - t τ 3 - 1.7796 * e - t τ
E F a s t t = 80.9888 + 150.7223 * e - t τ 1 + 140.2060 * e - t τ 2 + 94.0827 * e - t τ 3 - 0.1165 * e - t τ

The elastic solution and time-dependent elastic modulus values required for Alfrey’s correspondence principle are known.

4.1 Analysis 1 (Pugno-Carpinteri Method):

The geometric properties of the pipes used in modeling for this analysis, as well as the mechanical properties of both the pipes and the adhesive, are shown in Tables 1, 2 and 3, respectively.

Viscoelastic analysis proves to be more effective for flexible adhesives. It has been observed that joints bonded with flexible adhesive exhibit a more uniform stress distribution compared to those bonded with stiff adhesives. In Figure 12a and 12b, since E1 A1 < E2 A2, the maximum stress occurred at the end of the stiffeer material pipe (x= -c).

Figure 12
a, b. Shear stress distributions in Sika Power 4720 and Sika Fast 5215 adhesives over time were analyzed at intervals t=[0, 1, 10, 100, 1000, 10000] s, Mean shear stress: 6 Mpa

As the elastic modulus decreases over time, it was observed that the stress peaks migrated from the edges towards the middle of the overlap length. In Figure 12a, while the peak stress value decreased from 11.19 MPa to 9.89 MPa, the minimum stress value in the middle increased from 4.16 MPa to 4.61 MPa. Using the Pugno-Carpinteri method, changes in shear stresses along the midline of the adhesive were detected in the time interval t = [0 1 10 100 100 10000] s, and the results are shown in Figures 13a and 13b. Analysis has shown that the highest stress value occurs near the ends of the overlap length. The maximum stress peaks occurred at the edge of the stiffer pipe end. An observation was made indicating that at t = 0 s, there existed a notable dissimilarity between the shear stresses manifested within the central region of the connection and those experienced at its extremities. As time progressed, this dissonance in stress levels attenuated, leading to a more uniform distribution of shear stress across the joint.

Figure 13
a, b The variation of shear stress in the Sika Fast 5215 adhesive over time is examined at discrete time points t = [0, 1, 10, 100, 1000, 10000]. (Mean shear stress: 9,54 Mpa)

Subsequent to t = 10000 s, the shear stresses within the connection exhibited stability, with no further alterations observed. The maximum shear stress values of both adhesives are below the shear strength values.

Shear stress distribution at the shear strength limit of the adhesive:

At t=0 s, the shear strength value exceeded the safety limits, whereas at t=10000 s, it fell within the safety parameters. Notably, the shear stress value (at x = -10 mm and t = 0 s) surpassed the adhesive's shear strength value (11.19 > 10), rendering it the most critical site for cohesive crack initiation. Subsequently, the crack propagates over time, culminating in fracture. This finding aligns with (Abouel-Kasım, 2014).

Figure 14a, b and Table 4 illustrate the stress, strain and stiffness properties of the adhesive, organized according to shear modulus values at discrete time points of t = [0, 1, 10, 100, 1000, 10000] seconds, across two different configurtions. (See Appendix B Appendix B clc; clear all; subplot(1,2,1) % Figure 14a % Analytical Shear Stress Distribution under tensile load % Pugno and Carpinteri Method % Geometric and materials properties format long E1=207000; E2=207000;L=20;ta=0.2; F=12000; R2i=8 ; R2o=10; R1i=10.2; R1o=12.2; v1=0.3;v2=0.3;va=0.371; % Viscoelastic properties of adhesive Sika Fast 5215(Flexible) Ea1=151; Ea2=140; Ea3=94;Einf=81; eta1=104; eta2=1567; eta3=25190; E0=Ea1+ Ea2+ Ea3+ Einf; % Shear Stress Distributions A1= pi*(R1o^2-R1i^2); A2= pi*(R2o^2-R2i^2); beta=(E1*A1)/(E1*A1+E2*A2); t=[0.01 1 10 100 1000 10000]; n=length(t); figure(1) Smax=zeros(1,n); Gama=zeros(1,n); k=zeros(1,n); disp(' Table 4. Maximum shear stress and maximum shear strain variation '); fprintf(' t Ga Max. Shear Max. Strain Stiffness \n'); for i=1:1:n Ea(i)=E0*(1-((Ea1/E0)*(1-exp(-t(i)/(eta1/Ea1)))+ (Ea2/E0)*(1-exp(-t(i)/(eta2/Ea2)))+ (Ea3/E0)*(1-exp(-t(i)/(eta3/Ea3))))); Ga(i)=Ea(i)/(2*(1+va)); K=(2*pi*R2o*Ga(i))/ta; A1= pi*(R1o^2-R1i^2); A2= pi*(R2o^2-R2i^2); % alfa=(K*(E1*A1+E2*A2)/(E1*A1*E2*A2))^(1/2); alfa=sqrt(K*(E1*A1+E2*A2)/(E1*A1*E2*A2)); % beta=(E1*A1)/(E1*A1+E2*A2); beta=(E2*A2)/(E1*A1+E2*A2); c=L/2; A=exp(-alfa*c)/(exp(-2*alfa*c)- exp(2*alfa*c))+beta*((exp(alfa*c)-exp(-alfa*c))/(exp(-2*alfa*c)-exp(2*alfa*c))); B=exp(alfa*c)/(exp(2*alfa*c)-exp(-2*alfa*c))+beta*((exp(-alfa*c)-exp(alfa*c))/(exp(2*alfa*c)-exp(-2*alfa*c))); x=-10:0.01:10; % z=(x+c)/L; % E1*A1 < E2*A2 maximum shear stress close to x=-c shear=(-F/(2*pi* R2o))*alfa*(A*exp(alfa*x)-B*exp(-alfa*x)); plot(x,shear,'linewidth',2) title(' Shear Stress Distributions for Sika Fast 5215 NT', 'Color', 'b') xlabel('L- Overlap Length [mm]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k') ylabel('\tau_{rx} Shear Stress - [MPa]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); Smax(i)=max(shear); Gama(i)=Smax(i)/Ga(i); % The integral (area) of each curve obtained using the trapezoidal method gives the stiffness S(i)=trapz(x,shear); fprintf('%5.0f %3.f %3.4f %1.4f %3.4f\n', [t(i);Ga(i);Smax(i);Gama(i); S(i)]); hold on i=i+1; end text(-11.8,11.2,'t= 0 s','Fontsize',10, 'Fontweight', 'bold','Color', 'k'); text(-11.8,9.7,'t= 10000 s','Fontsize',10, 'Fontweight', 'bold','Color', 'k'); text(-6,11.5,'Areas under the curves are equal','Fontsize',12, 'Fontweight', 'bold','Color', 'r'); text(-4,11.3,[('K(0) = '), num2str(S(1),7)],'Fontsize',10, 'Fontweight', 'bold') text(-4,11.1,[('K(1) = '), num2str(S(2),7)],'Fontsize',10, 'Fontweight', 'bold') text(-4,10.9,[('K(10) = '), num2str(S(3),7)],'Fontsize',10, 'Fontweight', 'bold') text(-4,10.7,[('K(100) = '), num2str(S(4),7)],'Fontsize',10, 'Fontweight', 'bold') text(-4,10.5,[('K(1000) = '), num2str(S(5),7)],'Fontsize',10, 'Fontweight', 'bold') text(-4,10.3,[('K(10000) = '), num2str(S(6),7)],'Fontsize',10, 'Fontweight', 'bold') hold on x=linspace(-10,10,50); stress=F/(2*pi*R2o*L); stress_m= stress *ones(size(x)); Sm=trapz(x,stress_m); plot(x, stress_m,':b','linewidth',2) text(-4,9.7,[('Mean Stress \tau_m = '), num2str(max(stress_m),5)],'Fontsize',10,'Fontweight', 'bold','Color', 'b'); text(-4,8.8,[('K(\tau_m) = '), num2str(max(Sm),7)],'Fontsize',10, 'Fontweight', 'bold','Color', 'b'); hold on x=linspace(-10,10,50); strength=10; strength_m= strength *ones(size(x)); plot(x, strength_m,':r','linewidth',2) text(-4,10.1,'Shear Strength - 10 MPa ','Fontsize',10, 'Fontweight', 'bold','Color', 'r'); ylim([8.6 11.6]) text(-9.5,9,'K:Stiffness (N/mm) ','Fontsize',10, 'Fontweight', 'bold'); text(-8,8.8,'(a) ','Fontsize',10, 'Fontweight', 'bold'); grid minor hold on subplot(1,2,2) % Figure 14b yyaxis left semilogx(t,Smax,'-o','linewidth',2,'MarkerFaceColor','k') title('Shear Stress Distributions for Sika Power 4720'); xlabel(' t - [s]', 'Fontsize',12, 'Fontweight', 'bold'); ylabel('\tau_{rx} Maximum Shear Stress- [MPa]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k'); text(0.05,9.85,'(b) ','Fontsize',10, 'Fontweight', 'bold'); yyaxis right semilogx(t,Gama,'-o','linewidth',2,'MarkerFaceColor','k') title('Maximum Shear Stress / Strain Distributions ', 'Color', 'b'); xlabel('t - [s]', 'Fontsize',12, 'Fontweight', 'bold'); ylabel('\gamma_{rx} Maximum Shear Strain ', 'Fontsize',12, 'Fontweight', 'bold'); grid minor hold on semilogx(t,Gama,'linewidth',2) legend('Max.Shear Stress Sika Fast 5215','Max.Shear Strain Sika Fast 5215','Fontweight', 'bold','Location','Best'); Nomenclature a1, a, a2 Radii of the middle surface of the adhesive and Adherends A1, A2 Cross-sectional areas of adherends (pipes) D1, D2 Bending stiffness E1, E2 Young’s moduli of adherends Ea, Ga Young and shear moduli of the adhesive Ei Young’s modulus for the member ‘i’ ηi Viscous coefficient for the member ‘i’ F1, F2 Axial forces acting on pipes 1 and 2 at position x per unit length N1, N2 Hoop forces acting on pipes 1 and 2 at position x per unit length M1, M2 Longitudinal bending moments acting on pipes 1 and 2 at position-x Q1, Q2 Shear forces acting on pipes 1 and 2 at position x per unit length gi(z) Dimensionless functions gi, τi Prony series coefficients G∞ Equilibrium shear modulus of adhesive G(t) Shear relaxation module function H(t) Heaviside function (Unit step function) δ(t) Delta function G ̅(s) Shear relaxation module in Laplace domain KL Adhesive stiffness per unit lengthOverlap length (2c) P, Q Laplace Transform Polynomials R Outer radius of pipe 1(inner) r1i, r2o Outer radius of pipe 1(inner), inner radius of pipe 2(outer), ta Thickness of adhesive t1, t2 Thicknesses of adherends u1, u2 Axial displacements of pipes 1 and 2 v1, v2 Radial displacements of pipes 1 and 2 υa Poisson’s ratios of adhesive υ1, υ2 Thickness of adhesive α, β, λ Coefficients B1, B2 Coefficients C1, C2 Coefficients Λi, Bi, Ki, γij Coefficients τm Mean shear stress acting at mid-thickness of the adhesive layer τrx Shear stress evaluated at mid-thickness of the adhesive layer (Tensile) σrr Peel stress x, r, θ Cylindrical coordinates ).

Figure 14
a, b General stiffness in Sika Fast 5215 adhesive at discrete time points t=[0 1 10 100 1000 10000] s.
Table 4
Maximum shear stress and maximum shear strain variation in Sika Fast 5215 adhesive.

The integral of each graph corresponds to the axial force per unit length, denoting the stiffness of the adhesive.

This relationship is derived from Equation 18, F2π*a=τm*L=Constant [N/mm] .

The overall in-plane stiffness of the joint remains unchanged, independent of the modulus rating of its adhesive. These findings are in accordance with (Kumar & Pandey, 2010).

Maximum stresses diminish over time, reaching an equilibrium state at constant strain (stress relaxation),

Δ (%) Maximum Shear stress, 13.66% (Reduction in shear stress).

Maximum strains escalate over time, attaining an equilibrium state at constant stress (creep compliance).

Δ (%) Maximum Shear strain, 403.63% (Increase in shear strain).

4.2 Analysis 2 (Lubkin-Reissner Method):

The material utilized for the pipes comprises S235JR (steel), while Sika Power 4720 serves as the adhesive. The geometric properties of the pipes employed in the modeling for this investigation, as well as the mechanical characteristics of both the pipes and the adhesive, are delineated in Tables 1 and 2, respectively.

The mean shear stress corresponding to the applied tensile load of 7615.208 N is computed as 6 MPa utilizing Equation 18.

Employing the Lubkin-Reissner method and conducting 3D finite element analyses, the variations in shear and peel stresses along the midline of the adhesive within the time interval t = [0, 10000] s were ascertained, with the results depicted in Figure 15a and 15b.

Figure 15
a, b At time intervals t=[0, 10000] s, shear and peel stress distributions in Sika Power 4720 adhesive were investigated using the Lubkin-Reissner and 3D finite element methods.

Upon reaching equilibrium in adhesive stress relaxation, the shear and peel stresses at the tips witnessed a decrease. It was observed that the shear stress within the adhesive surpassed the peak value of the normal (peel) stress. Furthermore, it was noted that the degradation in the normal stress distribution occurred at a slower rate. This observation aligns with the earlier study by (Delale & Edoğan, 1981).

The viscoelastic analysis revealed that the stiffness of the adhesive joint system could increase over time, particularly in the middle section. It was observed that the peaks of the stress curves shifted from the ends towards the middle of the overlap length over time. This finding is in line with the observations of (Cheng, Özsoy, & Reddy, 2013).

5 CONCLUSION

This study investigated the viscoelastic stress distributions along the central portion of the adhesive in cylindrical lap joints under tensile loading, employing analytical and numerical methods based on Alfrey's correspondence principle. The obtained results were subsequently compared.

As the Alfrey correspondence principle primarily relies on the elastic solution of the problem, a thorough understanding of the elastic solution, both analytically and numerically, was necessary. To this end, the Pugno-Carpinteri method and the Lubkin-Reissner method were provided for analytical analysis, along with practical examples utilizing the Abaqus software in finite element analysis. Subsequently, transformed elastic modulus variables were determined using the generalized Wiechert model and the Laplace transform. Viscoelastic solutions for the relevant application were derived by substituting the “transformed elastic modulus variables” for the elastic modulus variables in the stress equations within the elastic solution.

The stiffness of the adhesive remains constant based on the stress distribution. With the decrease in the elastic modulus value over time, the stress peaks shifted from the edges to the middle of the joint lap length. It was observed that while the stress values decreased at the edges of the adhesive,increased in the middle.

Viscoelastic analysis proves to be more effective for flexible adhesives. It has been observed that joints bonded with flexible adhesive exhibit a more uniform stress distribution compared to those bonded with stiff adhesives. At the onset of loading, approximately at t = 0 s, maximum stress is observed near the edges of the adhesive, attributed to the elevated elastic modulus value (indicative of brittleness). When this value surpasses the shear strength of the adhesive, it initiates crack formation, ultimately resulting in damage within the adhesive.

One of this contribution is also the validation of the analytical approaches, as direct replacements for a more expensive numerical study. Thus, it is important that the differences of the results at the interface are explained, since in other regions the analytical models are validated. Additionally, Analytical solutions do not require special software or hardware so they can be used by anyone, anywhere. These advantages render analytical solutions superior to numerical and experimental methods.

Appendix A

clc;

clear all;

% Analytical Shear and Peel Stress Distributions under tensile load

% Geometric and materials properties

format short

E1=207000; E2=207000;Ea=2700;ta=0.2;L=20;F=1269.203;

R1i=8; R1o=10;R2i=10.2; R2o=12.2; v1=0.3;v2=0.3;va=0.367;

t1=2;t2=2;

Ga=Ea/(2*(1+va));

% Lubkin Reisner Method

syms g_1_0 g_1_0_t g_2_0 g_2_0_t g_3_0_t

syms t s

% Inner pipe properties (Steel).

L=20; c=L/2; a1=(R1i+R1o)/2 ;t1=2; E1=207000; v1=0.3;

% Outer pipe properties(Steel).

L=20; c=L/2; a2=(R2i+R2o)/2 ;t2=2; E2=207000; v2=0.3;

% Adhesive properties.

a=10.2;ta=0.2; Ea=2700; va=0.367; F=1269.203; tau=1;

Ga=Ea/(2*(1+va));

K1=(12*(1-v1^2)*(2*c/t1)^4*(t1/a1)^2)^(1/4);

K2=(12*(1-v2^2)*(2*c/t2)^4*(t2/a2)^2)^(1/4);

B1=((1-v1^2)*(2*c/t1)^2*((t1*Ga)/(ta*E1)))^(1/2);

B2=((1-v2^2)*(2*c/t2)^2*((t2*Ga)/(ta*E2)))^(1/2);

g11=(12*(1-v1^2)*(2*c/t1)^3*(2*c/t1)*((a*Ea*t1)/(a1*E1*ta)))^(1/4); % g11=gamma11;

g22=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t2)*((a*Ea*t2)/(a2*E2*ta)))^(1/4); % g22=gamma22;

g12=(12*(1-v2^2)* (2*c/t2)^3* (2*c/t1)*((a*Ea*t2)/(a1*E2*ta)))^(1/4); % g12=gamma12;

g21=(12*(1-v1^2)*(2*c/t1) ^3*(2*c/t2)*((a*Ea*t1)/(a2*E1*ta)))^(1/4); % g21=gamma21;

L1=2*v1*(2*c/t1)^2*(t1/a1); % L1=lambda1

L2=2*v2*(2*c/t2)^2*(t2/a2); % L2=lambda2

g30=-1; g20_2t=0; % fixed values at z=0

% T1=g10^”-3*a/a1*g30;

T1=sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2))-3*(a/a1)*(g30);

% T2=g10^”’-3*(a/a1)*g30^’

T2=sqrt(2)*K1*(sqrt(2)*K1*g_1_0_t-(K1^2)*g_1_0-2*sqrt(3)*v1*(a/a1)*(1/sqrt(1-v1^2)))-(K1^2)*g_1_0_t;

% T3=g20^”-3*(a/a2)*g30;

T3=g20_2t-3*(a/a2)*g30;

% T4=g20^”’-3*(a/a2)* g30^’=0;

T4=0;

% T5=g30-(B1^2)*g_1_0-(B2^2)*g_2_0;

T5= g30-(B1^2)*g_1_0-(B2^2)*g_2_0;

% T6=g30^’-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t;

T6=g_3_0_t-(B1^2)*g_1_0_t-(B2^2)*g_2_0_t;

A=[s^4+K1^4+g11^4,-g12^4,-3*(a/a1)*(s^2+L1);-g21^4,s^4+K2^4+g22^4,-3*(a/a2)*(s^2-L2);-(B1^2)*(s^2+L1), -(B2^2)*(s^2-L2),s^2-(a/a1)*(B1^2)-(a/a2)*B2^2];

B=(1/s)*[g_1_0*s^4+g_1_0_t*s^3+T1*s^2+T2*s-3*L1*(a/a1);...

g_2_0*s^4+g_2_0_t*s^3+T3*s^2+T4*s-3*L2*(a/a2);...

T5*s^2+T6*s-(a/a1)*(B1^2)+(a/a2)*B2^2];

G=A\B;

g=ilaplace(G);

disp(g)

g3z=g(3);

g1z=g(1);

g2z=g(2);

% Two function and three function derivatives unknown

% By arranging the g functions according to the 5 boundary conditions given in z=1, a 5x5 matris equation is formed. [M]*{D}=[P]

M=zeros(5);

q=subs(g3z,t,1);

for i=1:5

W=zeros(1,5);

W(i)=1;

M(1,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

end

P(1,1)=1-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

q=subs(diff(g1z,t,2),t,1);

q=subs(diff(g1z,t,2),t,1);

for i=1:5

W=zeros(1,5);

W(i)=1;

M(2,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

end

P(2,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

q=subs(diff(g1z,t,3),t,1)-3*(a/a1)*subs(diff(g3z,t),t,1);

for i=1:5

W=zeros(1,5);

W(i)=1;

M(3,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

end

P(3,1)=-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

q=subs(diff(g2z,t,3),t,1)+sqrt(2)*K2*subs(diff(g2z,t,2),t,1)+(K2^2)*subs(diff(g2z,t),t,1)-3*(a/a2)*subs(diff(g3z,t),t,1);

for i=1:5

W=zeros(1,5);

W(i)=1;

M(4,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

end

P(4,1)=0-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

q=subs(diff(g2z,t,2),t,1)+sqrt(2)*K2*subs(diff(g2z,t),t,1)+(K2^2)*subs(g2z,t,1);

for i=1:5

W=zeros(1,5);

W(i)=1;

M(5,i)=double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},W))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

end

P(5,1)=(-2*sqrt(3)*v2*((a/a2)*(1/sqrt(1-v2^2))))-double(subs(q,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},zeros(1,5)));

% g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t values

% The D values are substituted in the g functions and the actual g1z, g2z, g3z functions are found.

D=M\P;

z=linspace(0,1,20);

g=subs(g,{g_1_0, g_1_0_t, g_2_0, g_2_0_t, g_3_0_t},D'); % D' column vector

g1z=g(1);

g2z=g(2);

g3z=g(3);

figure(1)

% Plotting g1.

for i=1:length(z)

g1(i)=double(subs(g1z,t,z(i)));

end

plot(z, g1,'-k','linewidth',2)

title(' g_1(z) Dimensionless function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b')

xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

ylabel('g_1(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

legend('g_1(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','best');

text(0.1,-0.019,'(a) ','Fontsize',10, 'Fontweight', 'bold');

grid minor

figure(2)

% Plotting g2.

for i=1:length(z)

g2(i)=double(subs(g2z,t,z(i)));

end

plot(z, g2,'-r','linewidth',2)

title(' g_2(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b')

xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

ylabel('g_2(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

legend('g_2(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast');

grid minor

text(0.1,-0.018,'(b) ','Fontsize',10, 'Fontweight', 'bold');

figure(3)

% Plotting g3.

for i=1:length(z)

g3(i)=double(subs(g3z,t,z(i)));

end

plot(z, g3,'-k','linewidth',2)

title(' g_3(z) Dimensionless Function ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b')

xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

ylabel('g_3(z) Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

legend('g_3(z) Function','Location','Best', 'fontsize',10, 'Fontweight', 'bold','Location','NorthEast');

grid minor

text(0.1,-0.8,'(c) ','Fontsize',10, 'Fontweight', 'bold');

figure(4)

% Plotting Shear Stress.

for i=1:length(z)

g3(i)=double(subs(g3z,t,z(i)));

shear(i)=double(subs(diff(g3z,t),t,z(i)))/2;

end

plot(z, shear,'-r','linewidth',2)

title('\tau_{rx}- Shear Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b');

xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

ylabel('\tau_{rx}/ \tau_m Normalized Stress ', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k')

legend('\tau_{rx} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold');

xlim([0 1])

text(0.3,1.45,[('Lubkin/Reissner = '), num2str(max(shear),4)],'Fontsize',10, 'Fontweight', 'bold')

text(0.1,0.7,'(d) ','Fontsize',12, 'Fontweight', 'bold');

grid minor

figure(5)

% Plotting Peel Stress.

for i=1:length(z)

g2(i)=double(subs(g2z,t,z(i)));

v_2(i)=tau*((1-v2^2)/(E2*t2^2))*(2*c)^3*g2(i);

g1(i)=double(subs(g1z,t,z(i)));

v_1(i)=tau*((1-v1^2)/(E1*t1^2))*(2*c)^3*g1(i);

peel(i)=(Ea/ta)*(v_2(i)- v_1(i));

end

plot(z, peel,'-r','linewidth',2)

title('\sigma_{rr} - Peel Stress Function', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'b');

xlabel('z=x/L', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

ylabel('\sigma_{rr} / \tau_m Normalized Peel Stress', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k')

legend('\sigma_{rr} - Lubkin-Reissner','Location','northeast', 'fontsize',10, 'Fontweight', 'bold');

text(0.3,0.75,[('Lubkin/Reissner = '), num2str(max(peel),4)],'Fontsize',10, 'Fontweight', 'bold')

text(0.1,-0.3,'(e) ','Fontsize',12, 'Fontweight', 'bold');

grid minor

Appendix B

clc;

clear all;

subplot(1,2,1) % Figure 14a

% Analytical Shear Stress Distribution under tensile load

% Pugno and Carpinteri Method

% Geometric and materials properties

format long

E1=207000; E2=207000;L=20;ta=0.2; F=12000;

R2i=8 ; R2o=10; R1i=10.2; R1o=12.2; v1=0.3;v2=0.3;va=0.371;

% Viscoelastic properties of adhesive Sika Fast 5215(Flexible)

Ea1=151; Ea2=140; Ea3=94;Einf=81; eta1=104; eta2=1567; eta3=25190;

E0=Ea1+ Ea2+ Ea3+ Einf;

% Shear Stress Distributions

A1= pi*(R1o^2-R1i^2);

A2= pi*(R2o^2-R2i^2);

beta=(E1*A1)/(E1*A1+E2*A2);

t=[0.01 1 10 100 1000 10000];

n=length(t);

figure(1)

Smax=zeros(1,n);

Gama=zeros(1,n);

k=zeros(1,n);

disp(' Table 4. Maximum shear stress and maximum shear strain variation ');

fprintf(' t Ga Max. Shear Max. Strain Stiffness \n');

for i=1:1:n

Ea(i)=E0*(1-((Ea1/E0)*(1-exp(-t(i)/(eta1/Ea1)))+ (Ea2/E0)*(1-exp(-t(i)/(eta2/Ea2)))+ (Ea3/E0)*(1-exp(-t(i)/(eta3/Ea3)))));

Ga(i)=Ea(i)/(2*(1+va));

K=(2*pi*R2o*Ga(i))/ta;

A1= pi*(R1o^2-R1i^2);

A2= pi*(R2o^2-R2i^2);

% alfa=(K*(E1*A1+E2*A2)/(E1*A1*E2*A2))^(1/2);

alfa=sqrt(K*(E1*A1+E2*A2)/(E1*A1*E2*A2));

% beta=(E1*A1)/(E1*A1+E2*A2);

beta=(E2*A2)/(E1*A1+E2*A2);

c=L/2;

A=exp(-alfa*c)/(exp(-2*alfa*c)- exp(2*alfa*c))+beta*((exp(alfa*c)-exp(-alfa*c))/(exp(-2*alfa*c)-exp(2*alfa*c)));

B=exp(alfa*c)/(exp(2*alfa*c)-exp(-2*alfa*c))+beta*((exp(-alfa*c)-exp(alfa*c))/(exp(2*alfa*c)-exp(-2*alfa*c)));

x=-10:0.01:10;

% z=(x+c)/L;

% E1*A1 < E2*A2 maximum shear stress close to x=-c

shear=(-F/(2*pi* R2o))*alfa*(A*exp(alfa*x)-B*exp(-alfa*x));

plot(x,shear,'linewidth',2)

title(' Shear Stress Distributions for Sika Fast 5215 NT', 'Color', 'b')

xlabel('L- Overlap Length [mm]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k')

ylabel('\tau_{rx} Shear Stress - [MPa]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

Smax(i)=max(shear);

Gama(i)=Smax(i)/Ga(i);

% The integral (area) of each curve obtained using the trapezoidal method gives the stiffness

S(i)=trapz(x,shear);

fprintf('%5.0f %3.f %3.4f %1.4f %3.4f\n', [t(i);Ga(i);Smax(i);Gama(i); S(i)]);

hold on

i=i+1;

end

text(-11.8,11.2,'t= 0 s','Fontsize',10, 'Fontweight', 'bold','Color', 'k');

text(-11.8,9.7,'t= 10000 s','Fontsize',10, 'Fontweight', 'bold','Color', 'k');

text(-6,11.5,'Areas under the curves are equal','Fontsize',12, 'Fontweight', 'bold','Color', 'r');

text(-4,11.3,[('K(0) = '), num2str(S(1),7)],'Fontsize',10, 'Fontweight', 'bold')

text(-4,11.1,[('K(1) = '), num2str(S(2),7)],'Fontsize',10, 'Fontweight', 'bold')

text(-4,10.9,[('K(10) = '), num2str(S(3),7)],'Fontsize',10, 'Fontweight', 'bold')

text(-4,10.7,[('K(100) = '), num2str(S(4),7)],'Fontsize',10, 'Fontweight', 'bold')

text(-4,10.5,[('K(1000) = '), num2str(S(5),7)],'Fontsize',10, 'Fontweight', 'bold')

text(-4,10.3,[('K(10000) = '), num2str(S(6),7)],'Fontsize',10, 'Fontweight', 'bold')

hold on

x=linspace(-10,10,50);

stress=F/(2*pi*R2o*L);

stress_m= stress *ones(size(x));

Sm=trapz(x,stress_m);

plot(x, stress_m,':b','linewidth',2)

text(-4,9.7,[('Mean Stress \tau_m = '), num2str(max(stress_m),5)],'Fontsize',10,'Fontweight', 'bold','Color', 'b');

text(-4,8.8,[('K(\tau_m) = '), num2str(max(Sm),7)],'Fontsize',10, 'Fontweight', 'bold','Color', 'b');

hold on

x=linspace(-10,10,50);

strength=10;

strength_m= strength *ones(size(x));

plot(x, strength_m,':r','linewidth',2)

text(-4,10.1,'Shear Strength - 10 MPa ','Fontsize',10, 'Fontweight', 'bold','Color', 'r');

ylim([8.6 11.6])

text(-9.5,9,'K:Stiffness (N/mm) ','Fontsize',10, 'Fontweight', 'bold');

text(-8,8.8,'(a) ','Fontsize',10, 'Fontweight', 'bold');

grid minor

hold on

subplot(1,2,2) % Figure 14b

yyaxis left

semilogx(t,Smax,'-o','linewidth',2,'MarkerFaceColor','k')

title('Shear Stress Distributions for Sika Power 4720');

xlabel(' t - [s]', 'Fontsize',12, 'Fontweight', 'bold');

ylabel('\tau_{rx} Maximum Shear Stress- [MPa]', 'Fontsize',12, 'Fontweight', 'bold', 'Color', 'k');

text(0.05,9.85,'(b) ','Fontsize',10, 'Fontweight', 'bold');

yyaxis right

semilogx(t,Gama,'-o','linewidth',2,'MarkerFaceColor','k')

title('Maximum Shear Stress / Strain Distributions ', 'Color', 'b');

xlabel('t - [s]', 'Fontsize',12, 'Fontweight', 'bold');

ylabel('\gamma_{rx} Maximum Shear Strain ', 'Fontsize',12, 'Fontweight', 'bold');

grid minor

hold on

semilogx(t,Gama,'linewidth',2)

legend('Max.Shear Stress Sika Fast 5215','Max.Shear Strain Sika Fast 5215','Fontweight', 'bold','Location','Best');

Nomenclature
a1, a, a2 Radii of the middle surface of the adhesive and Adherends
A1, A2 Cross-sectional areas of adherends (pipes)
D1, D2 Bending stiffness
E1, E2 Young’s moduli of adherends
Ea, Ga Young and shear moduli of the adhesive
Ei Young’s modulus for the member ‘i’
ηi Viscous coefficient for the member ‘i’
F1, F2 Axial forces acting on pipes 1 and 2 at position x per unit length
N1, N2 Hoop forces acting on pipes 1 and 2 at position x per unit length
M1, M2 Longitudinal bending moments acting on pipes 1 and 2 at position-x
Q1, Q2 Shear forces acting on pipes 1 and 2 at position x per unit length
gi(z) Dimensionless functions
gi, τi Prony series coefficients
G Equilibrium shear modulus of adhesive
G(t) Shear relaxation module function
H(t) Heaviside function (Unit step function)
δ(t) Delta function
G ̅(s) Shear relaxation module in Laplace domain
K
L
Adhesive stiffness per unit length
Overlap length (2c)
P, Q Laplace Transform Polynomials
R Outer radius of pipe 1(inner)
r1i, r2o Outer radius of pipe 1(inner), inner radius of pipe 2(outer),
ta Thickness of adhesive
t1, t2 Thicknesses of adherends
u1, u2 Axial displacements of pipes 1 and 2
v1, v2 Radial displacements of pipes 1 and 2
υa Poisson’s ratios of adhesive
υ1, υ2 Thickness of adhesive
α, β, λ Coefficients
B1, B2 Coefficients
C1, C2 Coefficients
Λi, Bi, Ki, γij Coefficients
τm Mean shear stress acting at mid-thickness of the adhesive layer
τrx Shear stress evaluated at mid-thickness of the adhesive layer (Tensile)
σrr Peel stress
x, r, θ Cylindrical coordinates

References

  • Abouel-Kasem, A. (2014). Analysis and design of viscoelastic adhesively bonded tubular joint. Journal of Enginering and Apllied Sciences, Vol.1, 13-23.
  • Adams, R., & Peppiart, N. (1977). Stress analysis of adhesive bonded tubular lap joints. J. Adhesion, Vol.9, 1-18.
  • Aimmanee, S. (2018). A unified analysis of adhesive bonded cylindrical coupler joints. Applied adhesive bonding in science and technology, 11-31.
  • Alwar, R., & Nagaraja, Y. (1980). Viscoelastic analysis of an adhesive bonded lap joints. Computers & Structures, Vol.11, 621-627.
  • Arifoğlu, U. (2021). Diferansiyel Denklemler ve Matlab Çözümleri. Istanbul: Alfa BasımYayım.
  • Brinson, H., & Brinson, C. (2008). Polymer Engineering Science and Viscoelasticity An Introduction. New York: Springer Science-Business Media, LLC.
  • Campilho, R., Pinheiro, A., Moreira, R., & Sanchez-Arce, I. (2022). Validation of theoretical models for the strengh prediction of tubular adhesive joints. Procedia Structure Integrity, 41:60-71.
  • Cheng, F., Özsoy, Ö., & Reddy, N. (2013). Finite element modeling of viscoelastic behaviour and interface damage in adhesively bonded joints. Scriver Publishing LLC, 23-46.
  • Delale, F., & Edoğan, F. (1981). Viscoelastic analysis of adhesively bonded joints. Journal of Applied Mechanical, Vol. 48, 331-338.
  • Dragoni, E., & Goglio, L. (2013). Adhesive stresses in axially loaded tubular bonded joints-part-I. International journal of adhesion and adhesives, Vol.47, 35-45.
  • Gunay, D., Aydemir, A., & Özer, H. (1998). Stress analysis in cylindrical adhesive lap joint of pressure vessels. 8th international machine design and production conference semptember 9-11, Ankara, 339-349.
  • Karlsson, P. (2014). Determination of viscoelastic properties of adhesives. Linnaeus University, Sweden: Master’s Thesis in Structural Engineering.
  • Khan, M., & Kumar, S. (2017). Mitigating edge effects in adhesively bonded composite tubular lap joints under axial, pressure and torsional loading via stiffness-tailoring. 21st International conference on composite materials, 20-25th August. Xi'an.
  • Kumar, S., & Pandey, P. (2010). Behaviour of bi-adhesive joints. (24).
  • Lee, E. (1954). Stress analysisi in viskoelastic bodies. Brown University, 183-190.
  • Lubkin, J., & Reissner, E. (1956). Stress distribution and design data for adhesive lap joints between circular tubes. ASME, 1213-1221.
  • Machiraju, C., Phan, V., Pearsal, A., & Madanagopal, S. (2006). Viscoelastic studies of human subscapularis tendon: Relaxation test and a Wiechert model. Computer Methods and Programs in Biomedicine. Computer Methods and Programs in Biomedicine, 29–33.
  • Mathwork Inc. (2023). The student edition off Matlab.
  • Nagaraja, Y., & Alwar, R. (1976). Viscoelastic analysis of an adhesive tubular joint. Journal of Adhesion, Vol.8, 76-92.
  • Özer, H., & Çay, M. (2022). Viscoelastic model to evaluate the shear stresses in the bi-adhesively bonded lap joints. The Brazilian Society of Mechanical Sciences and Engineering, Vol. 44, 518.
  • Pugno, N., & Carpinteri, A. (2003). Tubular adhesive joints under axial load. Journal of Applied Mechanics, ASME, Vol. 70, 832-839.
  • Saraç, İ. (2021). Boru yapıştırma bağlantılarında farklı tasarım parametrelerinin yapıştırıcı tabakasında gerilme dağılımına etkisinin sayısal olarak araştırılması. Uludağ Üniversitesi Mühendislik Fakültesi Dergisi, 26(1).
  • Schapery, R. (1961). Two simple approximate methods of laplace transform inversion for viscoelastic stress analysis. Pasadena, California: Aeronautical Research Laboratory Wright Air Development Division.
  • Shishesaz, M., & Reza, A. (2013). The effect of viscoelasticity of polymeric adhesives on shear stress distribution in a single lap joint. The Journal of Adhesion, 89: 859-880.
  • Tapar, M. H., & Reza, A. (2021). Long- term shear stress distribution in adhesively bonded tubular joints under tensile load using the time- temperature superposition principle. The Journal of Adhesion, Vol. 97, 328-345.
  • Turner, A. J. (1948). Mechanical Behaviour of High Polymers. New York: Intescience Publishers.
  • Yilmazyurt M., Eyüpreisoğlu, S., Okyar, A. F., & Namlı, O. C. (2022). Viscoelastic characterization and mechanical hysteresis of commercial grade polypropylene. Journal of Polytechnic.

    Edited by

    Editor: Rogério José Marczak

    Publication Dates

    • Publication in this collection
      09 Aug 2024
    • Date of issue
      2024

    History

    • Received
      28 Apr 2024
    • Reviewed
      05 June 2024
    • Accepted
      05 June 2024
    • Published
      10 June 2024
    Individual owner www.lajss.org - São Paulo - SP - Brazil
    E-mail: lajsssecretary@gmsie.usp.br