Residual strength evaluation by DBEM for a cracked lap joint

The present work summarizes a numerical procedure aimed at the evaluation of the residual strength of a cracked lap joint, based on the competing failure mechanisms regulated by the R-curve analysis and plastic collapse. The model adopted for Stress Intensity Factors (SIFs) evaluation is based on the use of the Dual Boundary Element Method (DBEM) within the theoretical frame of Linear Elastic Fracture Mechanics (LEFM). The value of failure load was available from experiments, allowing a comparison with numerical results and consequent validation of the described procedure.


INTRODUCTION
he residual strength of an aircraft structure degrades during the life due to fatigue cracks, especially in the presence Multiple Site Damage (MSD), eventually evolving towards Widespread Fatigue Damage (WFD).Broek [1] published a report on the residual strength behaviour of 2024-T3 Alclad sheet panels.Small MSD cracks in combination with a lead crack can significantly reduce the load level driving to unstable crack propagation.Test data on residual strength of various types of airframes with multiple-site fatigue cracks are presented in [2], where it is showed how residual strength is affected by structural design features, bending stresses, material plasticity, arrangement of multiple-site cracks and stable growth of cracks under static loading.Structures with MSD fail when Stress Intensity Factors (SIFs) are within a range from the plane-strain fracture toughness (KIc) to plane-stress fracture toughness (Kc) and when net stresses are 30 to 90% of the yield stress.Full-scale tests show that the presence of MSD adjacent to a lead crack reduces the residual strength by 15% [3].SIFs can attain critical values in a way similar to strain energy release rates.The criterion for failure due to unstable crack growth can therefore be written as where KIc and Kc are the fracture toughness of the material under plane strain and plane stress conditions respectively.It is experimentally found that K Ic is constant for thick sections of a given material.K c is found to vary with crack length and component geometry and is applicable to thinner sections where stable crack growth can occur.From Eq.1, the failure criterion can be written as: T where Y is the geometric factor,  c is the critical stress and c a is the critical crack length at failure.This equation can be used to assess failure criterion for a component.Hence if a particular crack length is chosen and Y and K c are both known then Eq. 3 holds true: The critical stress  c must not be exceeded by the operating stress if failure of the cracked component is to be avoided.
In most cases, the critical stress will decrease as the crack length becomes longer and this must be considered in the long term assessment of working stresses.
If a stress level  c is chosen then the critical crack length is given by Eq. 4: The critical crack length must be much higher than the minimum detectable crack length a min so that the component can be inspected for crack growth at regular intervals.The above criterion does not take into account stable crack growth which can occur in thin sections of some materials.Under these condition the crack will only grow if the load is increasing whereas if the load is constant the crack will stop.
dK dK da da .R curves have been derived for many materials; more information on R curves and their use can be found in [4][5][6].

PROBLEM DESCRIPTION
he lap joint proposed for numerical residual strength assessment is represented in Fig. 1, with an MSD experimental scenario that is taken from [7].In Fig. 2 intermediate experimental cracked configurations, up to reaching 393000 fatigue cycles, are reported, whereas the modelled scenario (Fig. 3) makes reference to the experimental configuration after 399620 fatigue cycles, when the fatigue cycling was stopped and the experimental residual strength test started [7].The plates are made of Al 2024-T3 with Young's modulus E=72000 N/mm 2 and Poisson's ratio =0.33.The model adopted is based on the usage of Dual Boundary Element Method (DBEM) [8][9][10] under the hypothesis of Linear Elastic Fracture Mechanics (LEFM).In correspondence of 399620 fatigue cycles the residual ligament between hole N. 6 and 7 and between holes N. 8-9-10 had been cut by the advancing cracks and new cracks had nucleated from left hand side of hole N. 6 and from right hand side of hole N. 10; the same residual ligament failure happened between holes N. 15 and 16 with related nucleation of a new crack from right hand side of hole N. 16.The cracked plate of the lap joint has been modeled, using the commercial code BEASY, by a 2D single plate with constant traction applied on one side and constrains in y direction applied on the rivet holes in order to model pin actions, whereas no constraints are present in x-direction in order to allow transversal plate shrinkage (Fig. 3).With such constraints, the longitudinal plate compliance in the overlapping area (between the two rows of rivets) is neglected whilst it is underestimated in transversal x-direction, introducing an element of approximation.In the critical cracked area, the pin action modeling has been improved by explicitly inserting such pins in the holes (instead of just applying displacement constraints on the hole boundary) and moving the constraints on the pin centre.In particular, traction and displacements continuity conditions are imposed on 180 degrees of the pin-hole interface area (the supposed contact area after loading application), whilst the remaining part is disconnected by internal spring of negligible stiffness (Fig. 3).By means of a convergence study, it has been assessed that only in cracked holes and in holes nearby cracks it is necessary to explicitly model pins (rivets), whereas the remaining holes can be simply constrained against y-translation, as already mentioned.
The secondary bending (Fig. 4) is not considered, but the level of approximation introduced is judged acceptable for such kind of problems due to relevant plate flexural weakening caused by long cracks on the verge to become unstable (see [6] for a discussion on this topic).Gap elements have also been introduced to better tackle contact conditions [12] but the solution improvement has been judged quite negligible (less than 2% variation on SIFs), except in case of very short cracks initiated from the holes, more sensitive to pin-hole contact conditions.For this reason, and due to the computational effort of a non-linear analysis, they have no longer been used.The J-integral technique is adopted for SIF's evaluation, being more stable than Crack Opening Displacement method against crack mesh variations [13].
33 integration points are used along the J-integral path (Fig. 3c) whereas the increment of accuracy with 66 points turn out to be negligible.The mesh used for the lap-joint is based on about 326 quadratic elements: a p-convergence study has been realized showing that cubic elements provide an accuracy improvement of less than 2% and that 2 quadratic elements per 90 degrees are enough on the cracked hole, except when very short cracks are present (in such case 3 elements are recommendable, possibly with a scaling ratio).
After link-up of the cracks between holes and consequent development of the main crack there is no longer load transfer through the pins in the central part of such main crack even in the remote case they should not break when reached by the propagating crack.
Still remaining in the theoretical framework of linear elastic fracture mechanics, SIFs evaluation can be improved by empirically taking into account the elastic-plastic effects by the Irwin correction.Such correction is useful in a residual strength analysis and suggests to prolong the considered cracks of a virtual quantity calculated as a characteristic dimension of the plastic zone at the crack tip.Alternatively a fully elastic plastic nonlinear analysis could be attempted to get more accurate results [14].

SOLUTION PROCEDURE
wo approaches have been proposed for failure assessment:  Plastic collapse prediction, based on Von Mises stress exceeding 385 MPa, the average of yield ( y=330 N/mm 2 ) and rupture stress ( u =440 N/mm 2 ), in large zones of ligament;  R-curve analysis for stable and unstable crack growth assessment.With reference to the latter, it is well known that the failure criterion for plane strain structure is not valid for the case of thin metal sheet structure, because of extensive slow stable growth, under monotonic loading, prior to instability and catastrophic failure.Here rather than a single material parameter, a material curve (R-curve or KR-curve), representing an infinity of potential failure points (the crack length at instability is not known a priori), is necessary to make an accurate failure prediction.In this case two criteria must be satisfied to get an unstable crack growth: In the R-curve diagram there are two important points: 1.K o is the minimum SIF to start the crack propagation; 2. Kc is the critical stress intensity factor (instability point).Ko (the point of initial crack propagation) is independent from the specimen thickness and has a constant value equal to nearly 30 MPa m 1/2 for the considered material [6] so that the main crack will start propagating with a remote load equal to 63 MPa (as a matter of fact the first iteration in Table 1 starts considering 63 MPa).On the contrary Kc is strongly influenced from the specimen thickness: thinner specimens give higher Kc values and consequently exhibit slower stable crack growth.A sufficiently thick specimen will result in full plane strain and K c will then be equal to K Ic .
In order to obtain a crack driving energy (or force) curve an iterative process is needed, which is based on the following steps (Fig. 5):  the load is monotonically increased by small steps and for each of them a linear elastic analysis is performed by DBEM to calculate SIFs (when two consecutive configuration have nearly the same crack configuration it is possible to avoid the DBEM analysis, imposing a linear variation of SIFs vs. loads);  at each step cracks are prolonged by a length dai that is provided by the R-curve, as a function of the SIFs determined at the previous step; moreover, in order to provide the Irwin correction for SIFs evaluation, when the plastic effects become significant, cracks are prolonged by a virtual length r y =r p (Eq. 6);  for each crack tip, the G-curve (crack propagation driving force) is drawn and superimposed to the R-curve in order to find out the instability point, as resulting from the conditions in Eqs.5;  during the steady crack propagation some cracks will reach a link-up condition (Fig. 6) with other cracks or holes, when the plastic zone at the crack tip together with the plastic zone of the approaching crack or hole respectively, are covering the remaining ligament (Swift criterion).

NUMERICAL RESULTS
he first link-up (Table 1) is obtained with a load of 68 MPa, that is sufficient to create a plastic zone: T covering the ligament between the hole N. 5 and the crack tip N.1 (left side of hole N. 6 in Fig. 3); Keq is the equivalent SIF, calculated as in [9] but in this problem nearly coincident with the mode I SIF.When a propagating crack reaches a hole we have supposed the nucleation of a small, not detectable crack of 1 mm length on the opposite hole side (Fig. 7), as a consequence of local high stress gradients on an aged lap-joint (399620 fatigue load cycles).
Starting from the 8 th iteration, in order to improve accuracy of K values, the Irwin correction has been adopted.The second link-up (Table 1) between crack tip N. 2 (right side of hole N.10) and hole N. 11 (Fig. 3) is obtained with a load of 143 MPa, as suggested by Von Mises stresses exceeding 385 MPa in a large part of ligament (Figs.7b-c); again we have supposed the nucleation of a small, not detectable crack of 1 mm length on the right side of the hole N.11 (Fig. 7d).
Then, with the same load of 143 MPa and considering the main crack configuration obtained after the second link-up, the R-curve analysis provides an unstable crack growth (Fig. 8 and Table 2) of crack tip N. 1 up to hole N. 4 (Fig. 7d) and subsequent lap joint failure due to plastic collapse (the residual ligament cannot stand anymore the applied load and undergoes extensive yielding).Namely, with reference to the crack tip N. 1, the G-curve becomes tangent to the R-curve, reaching higher values and higher gradients (Fig. 8) and causing the aforementioned unstable growth.The R-curve (Fig. 9) used has the following equation [7]:  The G-curve, superimposed to the R-curve, is obtained in the following way:  crack tip N. 1 has been automatically propagated for a certain number of increments in order to get the SIFs for a variable crack length;  each SIF is corrected (with the Irwin criterion) by artificially modifying the correspondence between SIFs and related crack increment, in particular by backward shifting the crack length for each step of a quantity rp.From Figs. 10a-b it is evident that with a load of 143 MPa, before instability of crack tip N.1, Von Mises stresses are less than 385 MPa in most part of the ligament, in such a way that a failure based on plastic collapse is still premature.Consequently the real mechanism of lap joint failure is primarily related to fracture instability, responsible for the residual ligament reduction up to a condition in which the plastic collapse becomes effective.The numerical result of fracture instability at 143 MPa is close to the experimental collapse load equal to 139 MPa (specimen N. 5 in [7]) with consequent validation of the proposed procedure.

CONCLUSIONS
he procedure presented turn out to be a very effective way to model the assembly and exhibits a satisfactory agreement with experimental results, very attractive run times and an easy pre-processing phase (the mesh generation is very easy because based on monodimensional elements).The main advantages of the proposed DBEM two-dimensional approach to lap joint modelling are [see also [15][16]: each layer can be considered as an individual two-dimensional structure; individual layers can be explicitly modelled and connected with rivets (in case of need this provide a way to enhance the accuracy with respect to the simplified approach adopted in this work); the determination of the SIFs is straightforward with the J-integral technique; rivets can be modelled as separate DBEM zones, interacting with the main zone by gap elements (in case of nonlinear contact analysis) or by interface spring of negligible stiffness (to disconnect the pin-hole interface) or by simply enforcing continuity conditions.

Figure 3 :
Figure 3: DBEM model of lap-joint in the proposed initial cracked configuration, with highlight of: hole numbering on the cracked raw (a); main crack (b); hole constraints, rivets, J-path around the crack tip and "internal springs" (c).

TFigure 5 :
Figure 5: Flow chart of residual strength assessment procedure.

Figure
Figure 7a-d: Close up of the main crack after first link-up on the left side of main crack (a); internal point position, in the configuration preceding the second link-up (b) and related Von Mises stresses on the residual ligament of the main crack right side (c); close up of the main crack after second link-up on the right side (d).

Figure 8 :
Figure 8: R-G stability diagram with the G-curve corresponding to a load of 143 MPa.
in MPa*m 1/2 and da in mm; in this R-curve the Irwin plastic correction is included and, analogously, SIFs calculated by BEM are obtained with allowance for Irwin correction;

b
Figure 10a-b: Von Mises stresses (a) calculated with a load of 143MPa on the internal points (b) between crack tip N. 1 and hole N. 4 (see Fig.7d

Table 2 :
KR, KG data for the crack tip N. 1 (KG calculated with Irwin correction).