Bonded Joints Design Aided by Computational Tool

: In order to aid the design of bonded joints, a computational tool named System of Analysis for Joints (SAJ) was developed. The software can analyze single and double lap bonded joint with composite-­composite or metal-­composite mate-­ rials as adherent parts. Thus, SAJ can calculate the stress distribution, loads and displacements. Their results were compared (cid:87)(cid:82)(cid:3)(cid:191)(cid:81)(cid:76)(cid:87)(cid:72)(cid:3)(cid:72)(cid:79)(cid:72)(cid:80)(cid:72)(cid:81)(cid:87)(cid:3)(cid:86)(cid:82)(cid:73)(cid:87)(cid:90)(cid:68)(cid:85)(cid:72)(cid:3)(cid:11)(cid:36)(cid:37)(cid:36)(cid:52)(cid:56)(cid:54)


INTRODUCTION
In the last years, the usage of composite materials as a primary structural element has increased. Some new aircraft designs, such as Airbus A380 and Boeing 787, use composite materials even for primary structural elements, such as wing spars and fuselage skins, achieving lighter structures without reducing of airworthiness. One way to assemble those parts is to apply a bonded joint, which possesses some advantages, for doi: 10.5028/jatm.2012.04032912 example: better fatigue endurance;; joining dissimilar materials;; better insulation;; smooth surface and lower weight. However, there are also negative aspects, for example: no possibility to disassemble the joints;; peeling stress should be minimized and the preparation of the surfaces for bonding must be done carestress distribution in bonded joints and many design parameters Regarding this scenario, many researchers have conducted studies about bonded joints, trying to predict their mechaniexperimental tests or hybrid approaches, combining theoretijoints were carried out by Volkersen in 1938 (Mortensen, 1998). Volkersen modeled the adhesive layer as shear spring ural effects (Mortensen, 1998). The results obtained by the researcher showed that the load transfer at the adhesive was not uniform. Goland and Reissner (1944) improved Volkersen's model by considering the adherent like beams and the adhesive like springs in tension and shear loading. That model could analyze the distribution of the transversal loading induced by the secondary moment, which occurs in single lap joints. Hart-Smith (1970;; 1973), based on the last models, materials. Hart-Smith simulated the adhesive response, using an elastic-plastic material model. It is important to mention the main design parameters in the structural performance of the joint too. However, some of those models assumed that the laminates were isotropic or symmetric (only with tension ing apparent properties for the laminates. In 1975, Yuceoglu and Updike improved Hart-Smith's model by considering the transversal shear effects in the adherent. They showed that ity in shear loading, the normal strains can affect the stress distribution in the adhesive. Thomsen (1992) showed that an increase in the overlap length reduces the stress level in the adhesive layer. The researcher concluded that the application of adhesive layer with lower elastic shear and tensile moduli is better to use identical or nearly identical adherents in bonded joints, i.e., adherents with similar stiffness. Frostig et al. (1998) proposed an approach using high-order theory. In their analysis, the adhesive was modeled as an elastic continuum media (2D and 3D). So, the adhesive could transfer normal stresses, in-plane and transversal shear stresses. These results showed that there was a great variation of the normal stress in the adhesive close to the edges.
applying the Classical Laminate Theory (CLT). For the initial analyses, the adhesive layer was considered as homogeneous, isotropic and linearly elastic, but the adhesives used at aeronautical industry have non-linear behavior. Thus, for the follow-up analyses, the Springer and Ahn's model considered the elastic-plastic response for the adhesives (Springer and Ahn, 1998). Mortensen (1998), in his Ph.D. thesis, presented a development of a computational tool for analysis of bonded joints, showing the equations and hypothesis for various types of joints. The author showed the solving process of differential equations, using the multi-segment method of integration. Four years later, Ganesh and Choo (2002) showed the effect of spatial grading of adherent elastic modulus on the peak of stress and stress distribution in the single lap joint. This effect decreases the stress peak and provides a more uniform shear stress distribution in the joint. After that, Belhouari et al. (2004) showed a comparison between single and double lap joint using the advantages of using a symmetric composite patch in order to repair the damage parts;; also, that double patch had lower stress when compared with single patch repair. Regarding the manufacturing process parameters, Seong et al. (2008) showed that an increase of bonding pressure during the manufacturing process leads to higher strength bonded joints. The researchhigher strength bonded joints, and the major failure mode for single lap composite-aluminum joints is the delamination of the composite adherent. In the same year, based on experimental analyses, Kim et al overlap length and observed that the normal and shear stresses at the adhesive region reduce when the overlap increases.
However, it is important to mention that other researchers preferred to analyze bonded joints by the Finite Element Method, for example: Charalambides et al. (1998); ; Goyal et al. (2008). Besides, some researchers compared analytical methods to numerical approaches. Zarpelon (2008) carried out analyses to evaluate the mechanical behavior of single and stepped joints. Zarpelon concluded that Mortensen's model produces good results for symmetric laminates, but for laminates with higher asymmetry, the model should be improved. hypothesis, which is not adequate to simulate the laminates Ribeiro, M.L., Tita, V. with warping effects. Agnieszka (2009) showed a numerical method, regarding the sensitivity to hydrostatic stress for prediction of the delamination initiation. This method allows simulating the failure in the joint (at the overlap region) and in the adherent. Silva et al. (2009a;; 2009b) showed an excellent bibliographic review about models for bonded joints and performed a comparison between the most important models, showing the advantages and limitations for each one. More element method applied to simulate bonded joints.
In many studies commented earlier, the researches proved behavior of bonded joints due to the interaction phenomena (for example: mechanical, physical and chemical interactions) between parts (adherents) and adhesives. Besides, the the stacking sequence of the laminate can change the structural performance of the bonded joints. Therefore, it can be concluded that the design and analysis of bonded joints is very complex problem, which there is not still a closed solution. Thus, new proposals to solve this problem are required. In order to aid the design of bonded joints, a computational tool named System of Analysis for Joints (SAJ) was developed, which is capable of analyzing not only single lap joint, but also double lap joint. These joints could be made of composite-composite materials or dissimilar materials, i.e., metal-composite. For both types of joints, SAJ can calculate the stress and strain distributions, loads and displacements. In order to evaluate the limitations and advantages of the SAJ, some analyses were performed, using case studies. Moreover (ABAQUS TM (ESAComp TM bonded joint design parameters on its mechanical behavior was carried out (overlap length;; type of joint;; adhesive elastic modulus and adhesive layer thickness).

THE COMPUTATIONAL TOOL - SAJ
The computational tool was developed to calculate the joint loads, displacements and adhesive/adherents stresses with non-linear effects (Ribeiro et al., 2010). This computational tool named SAJ was programmed in Matlab TM language and can analyze single and double lap bonded joints. In fact, SAJ can predict the mechanical behavior of composite-composite and metal-composite bonded joints. Regarding composite adherents, the computational tool can calculate the stress and strain distributions in each layer of the laminate.
To predict the mechanical behavior of bonded joints, SAJ layer thickness in case of composite adherents, mechanical properties for adherents and adhesives, joint dimensions of adhesive and adherents, as well as loads and boundary conditions. Using the input data, SAJ builds a set of differential equations, which is solved by Matlab TM . After that, SAJ shows forces, displacements and adhesive stresses by graphics and tabular format.

Mathematical formulation
A set of differential equations of the multi-domain boundary value problem is implemented in SAJ. In order to obtain this set of differential equations, the problem domain (bonded joint) is partitioned in three regions: one part for adherents only, other part for the bonded region (overlap region) and the last part, again, only for adherents. Figure  1 shows these subdivisions, coordinate system, as well as an example of boundary conditions and loads for single lap joint (Fig.  1a) and double lap joint (Fig.  1b). For each region, free body diagrams of equilibrium in an Thus, the set of differential equilibrium equations for single and double lap joints can be written. Based on these equations and applying the hypothesis that all derivatives in y direction are equal to zero (cylindrical bending hypothesis), as well as considering plane stress state and Kirchhoff's kinematic relations, it is possible to obtain the set of differential equations shown by Fig.  2 For composite adherents, some terms of the equations are calculated using the CLT and the symbols a ij , b ij and d ij at Fig. 2 correspond to components of laminate compliance matrix. The membrane, coupling and bending-torsion compliances are represented by a ij , b ij and d ij , respectively. Also, at the Fig.  3, the sub index comma with "x" means partial derivate in "x". Ribeiro   Using the free-body diagrams (Fig.  2) and CLT, the set of equations for the subdivisions out of overlap region (single or double lap joint) is presented in Eq. 1. The upper letter i (i = 1, 2, 3) indicates the adherent. (1) Regarding the double lap joint, the same procedure described before (free-body diagrams and CLT) is used to obtain the set of differential equations for double lap joint adherent 1 (Fig.  2). Equation 2 shows the set of differential equations for this region. (2) For the other regions in the overlap area (double and single lap joint), the resulting set of differential equations is presented in Eq. 3. In this case, the upper letter i (i = 1, 2) indicates the adherent for single lap joint or i (i = 2, 3) indicates the adherent for double lap joint.
( x It is important to mention that the adhesive is simulated as springs under tension or compression combined to shear stresses as shown by the following equations: where t i is the thickness of the i adherent (i = 1, 2, 3), t a is the adhesive thickness, x is the rotation at the x axis, u 0 is the middle plane displacement in x direction and v 0 is the middle plane displacement in y direction.
Regarding the boundary conditions, in general they are provided as forces, moments or displacements in the edges of the problem domain which will be described in more detail in the following sections.
Finally, the differential equations for each subdivision are solved using Matlab TM , which can deal with multi-domain boundary value problem, also each subdivision are divided in n parts (mesh). It is important to mention that the ESAComp TM has basically the same formulation used in SAJ, but the solving process is different.

Numerical analyses procedure
The numerical analysis starts after SAJ reads input data from and adhesives mechanical properties, ply thickness and orientation (in case of composite materials), adhesive thickness and the dimensions, as well as loads and boundary conditions. Based on these data, SAJ calculates the stiffness and the compliance matrix. For composite parts (adherents), the CLT is applied.
Knowing the joint type, the boundary conditions and the compliance matrix calculated in a previous step, SAJ builds the set of differential equations presented earlier. This set of equations comprises on a boundary value problem, which is solved by "bvp4c" Matlab TM function. After that, SAJ shows the results by graphics and tabular format (Fig.  3).
For a better understanding of the differences between the numerical methods used, a short description of multiple point shooting method used by ESAComp TM and the Matlab TM funcassumed to be well-known by the readers.

Multiple point segment method (ESAComp TM )
Multiple point segment method is used for the boundary value problem with several initial conditions. In this method the domain is subdivided in n parts.
This method starts with an approximation for the equation derivative in one side of the domain (x = 0) regardless the value in the other side (x = 1) (Fig.  4). Moreover, this method uses the fourth order Range-Kutta (Butcher, 2003) to solve the set of differential equations 1 to 3.
With the initial shoot for the derivative in x = 0 the solution of the set of differential equations in x = 1 is compared with the prescribed value in this position. If the value of the tolerance), other approximation for the derivative in x = 0 is used. This procedure is repeated until the solution converges to the desired value.
The multiple point segment method requires low computational cost;; on the other hand, this method only works for simple problems (Saha and Banu, 2007). More details can be seen in Appendix I.

Matlab "bvp4c" function
The Matlab function "bvp4c" is the Simpson method to solve boundary value problem. Shampine et al. (2006) described how this method can solve the boundary value problem (see Appendix II for more detail).
The difference between the multipoint segment method, that is a shooting method, and the "bvp4c" function are that the solutions y(x) are approximated for the entire interval and the boundary conditions are considered every time during the solution. This method requires a discretization of the domain

Evaluation of the computational tool (SAJ)
commercial software ABAQUS TM was developed to compare (single and double lap joint) uses a second-order hexahedron element with 20 nodes (C3D20) for adherents, the adhesives were modeled with a second-order hexahedron element with 20 nodes (C3D20) too. ABAQUS TM constraint function "tie" is used to join the adhesive and adherents in the overlap region. This constraint function transfers all degrees of freedom element model for single lap bonded joint and Fig.  5b shows Ribeiro, M.L., Tita, V. composite single and double lap bonded joints were investigated. The adherent and adhesive mechanical properties and the joint characteristics are shown in Table  1. A normal load in x direction of 1 N/mm was applied on single and double lap joint. This load is small enough in order to avoid inelastic strains in adherents or in adhesive. Also, all the stresses in the adhesive layer will be divided by the applied stress (normalized stresses) in order to improve the comprehension of the load transfer by the bonded joint.
z direction for single lap joint. It is possible to verify the difference between element model, at the edge of adherent 1, all displacements (x, y and z the opposite edge, at adherent 2, z Bonded Joints Design Aided by Computational Tool  rotations are free and the loading is applied in x direction (Fig.  5a). For SAJ, at the edge of the adherent 1, x(u) and z(w) y(v) is free (Fig. 1a). In rotation at y (v) is free and the normal loading is applied in x(u) direction (Fig.  1a).
and SAJ are almost the same. These small differences are due to limitations of the hypothesis adopted for SAJ and the mathematical procedure used to solve the set of differential is a little bit higher than SAJ, i.e. the SAJ model is slightly TM do not regard for free edge effects on the adhesive layer.
It is important to mention that the commercial software ESAComp TM has some limitations of boundary conditions once only three types of boundary conditions are available. In SAJ uses the same set of equations used in ESAComp TM , boundary conditions, but still has some limitations (for some boundary conditions was not possible to solve the set of differential equations). Despite the limitations to apply boundary conditions showed by ESAComp TM element model is possible to simulate any boundary conditions. So, in order to proceed with the analysis, the boundary conditions between SAJ and ESAComp TM were keep as Due to ESAComp TM limitations, the boundary conditions used for simulations were similar to SS.
The differences between SAJ and ESAComp TM responses (Fig.  6) are due to numerical differences to solve the set of differential equations and small differences in the boundary entire single lap joint. The joint dimensions and coordinates were shown in Fig. 1a for single lap joint and Fig. 1b for double lap joint. model, in the edge of adherents 2 and 3, all displacements and only the loading is applied in x direction (Fig.  5b). For SAJ, in the edge of the adherents 2 and 3, x(u) and z(w) displacements edge, at adherent 1, only the normal loading is applied in x(u) direction (Fig. 1b). The same consideration for boundary condi-Regarding to ESAComp TM , the CF boundary conditions were adopted which are, as for the last case, the most similar boundary conditions available in the software to proceed with the evaluations. In this case, despite the differences, the results are very close to SAJ (Fig.  6). Figure  6c shows the normalized az and ax for single lap bonded joint and all results are again similar too. In this case, it is possible to observe that the normal stress in the adhesive layer az has a relatively high tensile stress in the adhesive edge, than it changes to compressive stress in less than 2.5 mm and in 5.0 mm the normal stress is almost zero or has a very low value. The same behavior is observed in the other edge. Equation 6 explains this behavior (for SAJ and ESAComp TM ) once it regards the relative displacement has the same behavior for normal stress in the adhesive layer. Figure 6d shows that SAJ results of az and ax for double lap the ESAComp TM results are different. These differences depend on the boundary conditions imposed by each software, as well as the solution algorithm used in each model. It is important to note that the highest differences are observed at the end of the overlap region, where the highest stress values are observed. The stress singularity is not observed for double lap joint. For the second set of analyses, metal-composite single and double lap bonded joint were investigated. Again a normal load in x direction of 1 N/mm was applied on single and double lap joint. There is no change in the composite lay-up.
For single lap joint, aluminum was used for adherent 1 and composite laminate for adherent 2 (Fig.  1a). For double lap joint, aluminum was used for adherent 1 and composite laminates for adherents 2 and 3 (Fig.  1b). The material properties are shown in Table  1. The boundary conditions and loading for single and double lap joints were applied as commented earlier.
single lap joint, where ESAComp TM in the composite side and stiffer in the aluminum side. SAJ results are between ESAComp TM results. This occurs due to differences in computational method applied to solve the problem, as well as the boundary for the entire double lap joint and the SAJ results converge to TM results are different due to the applied boundary conditions. Figure  7c shows az and ax for single lap bonded joint and all results are very close. The stress singularity appears again for single lap hybrid joint.
For double lap joint, there are differences between the results of az and ax (Fig. 7d). These differences depend on the applied boundary conditions, as well as the solution algorithm used in each model, as discussed earlier. It is important to note, again, that the highest differences are observed at the end of the overlap region, where the highest stress values are TM are very similar, This occurs because SAJ and ESAComp TM use a cylindrical As commented before, the stress singularity that appears for single lap joints (composite-composite or composite-metal) is due to the relative displacements between the adherents (Eq.6). Once the set of equations (Eqs. 2 and 3) regards the adhesive stresses they affect the solution and then in the adherents kinematics. In the single lap joint solution, inside the overlap region, the w displacement of adherent 2 (see Eq. 6) becomes bigger (in absolute value) than the vertical displacement of adherent 1. This behavior leads to a compressive normal stress in the adhesive layer. Also it is important to mention that the numerical procedures used in SAJ or in ESAComp TM smooth the results between each subdivision avoiding discontinuthe relative displacement between the adherents leads to a compressive normal stress close to the adhesive layer edges.
Considering the results shown earlier, it is possible to conclude that SAJ is a computational tool, which can predict the mechanical behavior of composite-composite and metalcomposite bonded joints (single and double) with similar accuracy shown by other software. Therefore, in the next section, some design parameters will be investigated in order to joints.

Design parameters study for composite bonded joints using SAJ
One of the most important components of the joint is the adhesive layer at the overlap region, where physical and chemical interactions between adherents and adhesive occur and the load is transferred from the part (adherent) to the adhesive and vice versa. Thus, design parameters such as overlap length, type of joint, adhesive elastic modulus, and adhesive layer thickness, which affect the stress distribution in the overlap region, mainly the az and ax component stresses, were investigated. For all parametric studies, the boundary conditions and loading were applied as described here in section Evaluation of the computational tool (SAJ). The mechanical properties and other important characteristics for adhesive and adherents are given in Table  2 for both types of joints (single and double).

Effect of the overlap length
For this investigation, metal-composite single and double lap joints were used. A load of 1 N/mm was used for both types of joints. Aluminum 2024-T3 was assumed in single lap . For the double lap joint, aluminum 2024-T3 was assumed for adherent 1, and composite material was applied for other adherents. An epoxy adhesive was considered for both types of joints. Figure 8a shows the effect of the overlap length in the adhesive layer stress distribution for a single lap joint and Figure 8b shows the results for double lap joint. Table 2 shows the az and ax values obtained at the left edge of the adhesive layer for single and double lap joint (Fig.  1). Based on the results, it is observed that an increase in the overlap length leads decreases of stress state in the adhesive edge, mainly Ribeiro, M.L., Tita, V.  for single lap joints, considering the lengths studied in this paper. Besides, the rate for stress reduction decreases with an increase in the overlap length. So, it is possible to conclude that there is a length in which any further increases in the overlap length do not decrease the stress state in the adhesive layer edge. This trend is clearer for the double lap joint.

Comparison between single and double lap joints
For this study, all the joint parameters for single and double lap joint were kept the same as for the overlap effect analysis. It is important to note that the overlap length is equal 20 mm. The results are shown for half of the overlap length in the region with greater differences between these two types of joints. Figure  8 shows the difference between double and single lap joint stress distribution in the adhesive layer for the same load conditions and same joint characteristics. Figure 9a shows the difference between these two types of bonded joints for shear stress in the zx-plane ( ax ) and Fig.  9b shows the difference for the normal stress ( az ). As expected, the amplitude range of az and ax lap joints.

Effect of the adhesive elastic modulus
Another important parameter is the adhesive elastic modulus. This parameter was investigated, keeping other parameters constant and using three realistic values for the adhesive elastic modulus (1.5GPa;; 2.0GPa;; 3.0GPa) found in the literature (San Román, 2005).
Analysis of the effect of this parameter on the stress distribution in the adhesive layer was performed. In this case study, a normal load of 1 N/mm was imposed for both types of joints. Figure  10a shows the results for single lap joint, and Fig.  10b for double lap joint. It can be observed that adhesives with lower values of elastic modulus lead to lower stress state in the overlap region for single and for double lap bonded joints. This can be explained low strength values, and in a real joint design, it is desirable that the adhesive has a satisfactory performance. Thus, it is important to balance the stiffness and the strength of the adhesive layer. Besides, the adhesive can show inelastic strains according to the level of loading and the yielding stress of the polymer.

Effect of the adhesive layer thickness
The adhesive layer thickness affects the stress distribution in the adhesive layer. This parameter was investigated, keeping the other parameters constant, regarding three realistic adhesive layer thicknesses (0.05 mm;; 0.5 mm;; 1.0 mm), which are found in the literature (Qian and Sun, 2009). The mechanical properties and composite adherents and adhesive characteristics are shown in Table 1. The analyses were carried out using a normal load of 1 N/mm. This load leads to adhesive stresses low enough to avoid adhesive non-linear behavior, or any adhesive or adherent failure. Figure 11a shows a single lap joint adhesive stresses distribution for three different thicknesses, and Fig.  11b for double lap joint. These results show that this parameter can mainly close to the edge of the overlap region. Adhesives with lower values of thickness lead to higher stress state in the overlap region for single and for double lap bonded joints. This can be explained by Eqs. 4, 5 and 6.

CONCLUSIONS
SAJ, a new computational tool, has shown to be adequate in performing composite-composite and metal-composite single and double bonded joint analysis. Therefore, very quickly, it is possible to analyze a set of different joints, varying many parameters, for example: materials (adhesives and/or adherthickness of the laminate and/or the adhesive;; overlap length and the type of the joint (single or double). However, due to provided some small deviations in the results, when compared the analyses from SAJ are more conservative, which is very interesting for conceptual and preliminary design of a product.
Regarding to parametric study, SAJ leads to some conclusions, which can be used as a guide during product design. For example, a thicker adhesive layer (keeping other parameters constant) could reduce the adhesive edge stress state, increasing the strength of the joint, but thicker adhesives could lead to adhesive cohesive failure. The adhesive stiffness affects peak at the edges of the adhesive layer. Therefore, it is more it is important to verify what is more important to the product in service. Short overlap length increases the stress peak at the adhesive layer edges. Thus, it is reasonable to increase the overlap length, but the joint weight could increase too. Therefore, it is very important to balance all parameter values. Finally, by understanding the behavior of the stress distribution in the joint, it is possible to design products made of bonded joints with more accuracy even during the conceptual and preliminary phases. Ribeiro, M.L., Tita, V.

APPENDIX I
For Multiple Point Segment Method, consider the linear differential equations system in matrix notation The boundary conditions could be described as: The solution can be regarded as: Where: G is the integration constant, Y(x) is the general solution of the system and Z(x) is the particular solution. Consider: Solving equation 14 in x=a, yield G=Y(a) and for x=b . Applying this results in (10), yields Saha and Banu (2007) showed other example of the multiple point segment method application.

APPENDIX II
For the Matlab function "bvp4c", consider a differential equation y"+ y = 0 smooth in the [a,b]. The boundary condition in x = a is y(a) = A and y'(a) = s, the solution of the y(b,s) = . Consider that the algebraic equation s has a solution. Regard a function u(x) as the solution for y(a) = A and y'(a) = 0 and v(x) be the solution for y(a) = 0 and y'(a) = 1, this linear approach yields in y(x,s) = u(x) + sv(x) and with the boundary condition = y(b,s) = u(b) + sv(b), its results in a linear set of algebraic linear equations with initial derivative equal s.
In some kind of problems, the equations to solve the problem are non-linear. It implies that the existence and the In practice, to solve a set of differential non-linear equations are based in the solution of initial value problem and in non-linear equations solvers. Matlab TM "bvp4c" function uses the colocation method to solve the problem. For example, consider the following differential equation: ' , , , And the non-linear boundary conditions in x = a and x = b are: Where p is the vector of unknown parameters. The approximate solution, S(x), is a third order polynomial function smooth in the interval [x n , x n+1 ] with mesh a = x 0 < x 1 < … < x n < b, satisfying the boundary condition g(S(a), S(b)) = 0 the interval limits as well as the differential equations inside the mesh. This results in a non-linear system for S(x). Ribeiro, M.L., Tita, V.