Plasticity-induced crack closure identification during fatigue crack growth in AA2024-T3 by using high-resolution digital image correlation
Florian Paysan David Melching Eric Breibarth

Interacting Multiple Model-based Joint Homography Matrix and Multiple Object State Estimation
    [
    
=============================================================================================















§ ABSTRACT 
Fatigue crack growth in ductile materials is primarily driven by the interaction between damaging and shielding mechanisms. In the Paris regime, the predominant mechanism for retardation is plasticity-induced crack closure (PICC). However, some of the mechanisms behind this phenomenon are still unclear. Identifying and separating the three-dimensional aspect from other shielding aspects during experiments is extremely complex. In this paper, we analyze the crack opening kinematics based on local crack opening displacement measurements in both 2D  high-resolution digital image correlation data and 3D finite element simulations. The results confirm that the crack opening stress intensity factor K_op differs along the crack path. We present a new method to determine K_op at the crack front allowing us to identify PICC as the predominant shielding mechanism in fatigue crack growth experiments. Furthermore, this work contributes to the discussion on the damage-reducing effect of PICC, since we find that the influence on fatigue damage in the plastic zone remains negligible when the crack is closed and crack surface contact is directed towards the surface.  







§ INTRODUCTION
 


From a material science perspective fatigue crack growth is mainly driven by an interaction between damaging and shielding mechanism <cit.>. While damage mechanisms occur in front of the crack tip and lead to crack propagation, shielding mechanisms have a retarding effect. Elber <cit.> introduced the concept of crack closure in 1970. Under cyclic loading, the crack faces remain partially closed even though the fatigue crack is not completely unloaded and should be open. Assuming that no damage is induced in the process zone in front of the crack tip, he presents the Δ K_eff concept:

    Δ K_eff = K_max-K_op

K_op is declared as crack opening stress intensity factor and describes the crack tip stress, at which the crack is fully opened for the first time during loading. The crack closure is counted among the shielding mechanisms due to its damage-reducing effect according to Ritchie <cit.>. Elber <cit.> related crack surface contact to plastic deformation within the plastic wake of a growing fatigue crack. This phenomenon is also known as plasticity-induced crack closure (PICC) and is assumed to be the predominant crack closure mechanism in the Paris regime of ductile materials <cit.>. Furthermore, it is the only crack closure mechanism that can be fully explained physically <cit.>. For generic classification, the literature typically distinguishes between plane stress and plane strain conditions. Dugdale <cit.> and Newmann <cit.> explain PICC for plane stress conditions by necking near the crack tip during crack opening. 
The material deforms plastically towards the crack surface, which is also known as out-of-plane flow. As a result, contact between the crack surfaces occurs before the crack is fully unloaded. While this explanation is widely accepted in the fracture mechanics community <cit.>, the existence of PICC under plane strain conditions is debatable since, per definition, no material flow is allowed in the thickness direction <cit.>. Using finite element (FE) simulations, Riemelmoser and Pippan <cit.> recognize that material flows towards the crack tip while the crack opens. In contrast to PICC under plane stress conditions, this material flow only leads to contact right behind the crack front. 
Numerical investigations by other authors validate that PICC is mostly present in plane stress conditions near the free specimen surface <cit.>. 
However, 3D FE simulations are essential to investigate crack closure, but are very time-consuming due to the mesh refinement. Paysan  performed a parameter study as a function of the load and geometry to investigate the crack contact kinematics <cit.>. They conclude that the maximum stress intensity factor K_max and the sheet thickness t are the main factors that determine the location of the crack surface contact. However, there is still a need for further research, e.g. on how 3D crack surface contact affects damage development within the plastic zone or how other mechanisms such as crack deflections or branching interact with PICC.

The Δ K_eff concept requires the accurate determination of  K_op. The ASTM E647-15 standard <cit.> recommends the compliance method, which assumes that K is proportional to the crack opening <cit.>. Backface strain or crack mouth opening displacement (CMOD) gauges are commonly used for measurement (see Figure <ref>a). If crack closure is present, it results in a non-linear opening behavior. Both domains can be separated within the crack opening curve (see Figure <ref>b) enabling the determination of K_op.



Both the large distance to the crack tip position and the strong scattering in strain gauges are the main disadvantages of using compliance gauges <cit.>. Tong <cit.> claims that the study of crack closure from a systematic perspective became possible only with high-resolution digital image correlation (HR-DIC). Sutton <cit.> investigated the crack opening displacement (COD) in HR-DIC full-field data as early as 1983. Carroll <cit.> evaluate K_op based on COD measurements at different positions along the crack path and concluded that K_op strongly depends on the measurement position within the full-field displacement DIC data. Comparisons with FE simulations <cit.> confirm these findings and support the conclusions made by Carroll <cit.>. According to O'Conner <cit.> and Tong <cit.>, the value of K_op increases as the distance to the crack tip decreases. However, a systematic relation between K_op and measurement location has not yet been found. The non-linearity of the crack opening curve provides information about the crack contact <cit.>. If the transition to the linear domain is smoothly curved, a continuous opening crack is assumed. If the transition is abrupt, it indicates a sudden loosening of the crack surfaces. This behavior is often associated with roughness-induced crack closure. The impact of PICC on damage evolution or crack driving force has been a topic of debate for some time <cit.>. Research by Tong <cit.> questions Elber's <cit.> claim that no damage is caused when a crack is partially closed. It was found that the J-integral result is only slightly different when crack closure is not taken into account. Vasco-Olmo <cit.> followed their conclusions and proposed ΔCTOD_p as a damage describing parameter. They found that crack tip opening displacement (CTOD) measurements are characterized by another non-linear domain in the upper region of the crack opening curve (see Figure <ref>b). They correlate this behavior to the induction of plastic strain as a cause of the plastic zone evolution. Based on their numerical analysis, Oplt <cit.> correlated the ΔCTOD_p to the size of the plastic zone along the crack front. They concluded that ΔCTOD_p represents the plastic strain accumulation within the plastic zone region. A comparable methodology is employed by the δ_5 concept, which is predominantly utilized for the investigation of fatigue cracks in thin structures <cit.>. A clip gauge is employed for the measurement of displacements in close proximity to the crack tip, with the reference points at a distance of 5mm from each other.

However, the position the two strain gauges can remain unchanged with the advancing fatigue crack and its crack tip position. Consequently, the location dependency effect of the measurement position is not incorporated.

In the present work, we address the dependence of K_op on the measurement position, the influence of the 3D crack-surface contact on the damage evolution within the plastic zone, and the difficulty in identifying the dominant crack-closure mechanism in experimental data. We aim to improve the general understanding of PICC. Therefore, we use 3D FE crack propagation simulations and link them to crack opening curves based on HR-DIC data. The excellent agreement allows the investigation of the 3D contact characteristics and the associated damage accumulation within the plastic zone.




§ METHODOLOGY
 



 §.§ Specimen, Material and Loads

We conducted fatigue crack growth experiments according to the ASTM 647-15 standard <cit.>. A middle tension (MT) specimen with a width of W=160mm is used as a basis for investigating  fatigue crack growth. The dimensions of the specimen are given in Figure <ref>.



In this study, we consider the aluminum alloy AA2024-T3, a material with a wide range of applications in the field of aeronautics, as a reference material.
The specimen was loaded in the rolling direction, resulting in the fatigue crack growing perpendicular to the maximum grain elongation of the sheet material. The mechanical properties are given in the Table <ref> and taken from the literature <cit.>. 




The fatigue crack starts from an initial notch with a length of a_i=8mm. We investigate mode I fatigue crack propagation at sinusoidal load with a constant amplitude in y-direction. The maximum force F_max=15kN and the load ratio R=0.1 are constant. At a final crack length of a_f=27.8mm we collected data to investigate fatigue crack closure, where stress intensity factors of K_max=14.9MPa√(m) and Δ K=13.4MPa√(m) are present.  



 §.§ Finite Element Simulation

In general, the numerical simulations refer to the FE model presented in <cit.> with some crucial adaptions concerning loading constraints and meshing. All FE simulations are performed with ANSYS Mechanical APDL on a RedHat Linux workstation with two Intel Xeon Gold 6240 18C CPUs and 256GB memory (DDR4-2933 RAM). For the sake of comprehensibility, a brief summary of the FE model is following given:

The FE model is based on the geometry of the MT(160) specimen, shown in Figure <ref>. To reduce model complexity and computational resources, we exploit symmetries of the specimen to reduce the geometry to a 1/8 model. The clamping region is excluded and all holes are neglected enabling a structured mapped meshing. Figure <ref> shows the constraints definition. The load F is linked into a pilot node that is coupled to all nodes on the top of the model. In addition, the degrees of freedom of all nodes in the symmetry planes are restricted vertically to their corresponding plane. Within the x-z plane at the bottom of the model, the initial center crack with an initial crack length of a_i is defined. All symmetry constraints within the interval x ∈ [0, a_i] are deleted to allow free deformation of the crack surfaces. 
The following Figure <ref> sums up the most important features of the FE model.



In general, the model is meshed using a mapped mesh strategy with 8 nodes of linear hexagonal SOLID185 elements.
The structural meshing can be divided into three separate sections: (a) the global mesh of the overall model, (b) a more refined section, schematically shown in Figure <ref>a, in which the following analyses take place, and a transition section, which connects both. Within the refined section, the volume elements have an aspect ratio of 1 and an element size of e_FE,f=0.08mm. Especially contour and deformation analyses benefit from the chosen strategy, as the element shape in the observation space is always identical. However, the element size is much larger than recommended in literature <cit.>. Choosing smaller element sizes results in strong deformations of single elements, which leads to convergence issues in combination with the element contact definition applied. The mesh outside the refined section is characterized by a flexible element size of e_FE,g=1mm. The contact definition is realized via contact elements. The contact is assumed to be asymmetric and all surface effects, such as crack surface roughness, are neglected. As a counterpart to the freely deformable MT(160) model, a rigid target surface in the form of a plate (SOLID185, element size e_FE,c=1mm) is defined. Following, a pairwise frictionless rigid-to-flexible contact formulation is set between the two model components. The solid elements at the crack surface of the MT(160) model are covered with CONTA174 elements and TARGE170 elements overlay the elements of the rigid contact body.  Furthermore, we use the Augmented Lagrange Algorithm as a contact solver. The contact stiffness and penetration are set to k_c=1N/mm^2 and z_c=1μ m with the aim of reducing the contact penetration as much as possible. The contact definition is completed with the instruction that all initial gaps are closed. The implemented model reflects the elasto-plastic material behavior using bilinear isotropic hardening and is based on the mechanical properties of AA2024-T3 given in Table <ref>. Anisotropic properties are neglected within the FE simulation. 

In order to analyze PICC, crack propagation needs to be performed within the elastic-plastic FE simulation. This causes a plastic wake leading to crack surface contact. We use the Releasing-Constraint <cit.> method to conduct crack propagation. A single load cycle consists of a loading and unloading step followed by a node-releasing step, in which all constraints of one element row in z-direction are  released. Consequently, the crack advances exactly one element length e_FE,f within the refined mesh with each cycle at minimum load F_min. Taking N=50 cycles into account, a total crack propagation of Δ a= a_f-a_i=4mm is achieved. The application of only one loading and unloading step between crack propagation is in contrast to the general recommendation in the literature <cit.> to iterate this at least twice. The authors <cit.> assume that the influence of cyclic plasticity can be neglected and does not have a significant impact on the values of K_op as reported by Camas <cit.>. 

When the total crack advancement Δ a=4mm and cycle number N=50 is reached, another load step is added, which is subdivided into 100 sub steps allowing for a detailed characterization of the crack opening behavior.



 §.§ Fatigue crack growth experiments

In order to compare the numerical results with experiments, we perform fatigue crack growth in MT(160) specimen made of AA2024-T3. The experimental setup is illustrated in Figure <ref>, with a servo-hydraulic testing machine inducing the fatigue loading in terms of sinusoidal load with constant amplitude. Load parameters are equal to the numerical studies F_max=15kN at R=0.1. For analyzing crack closure characteristics, we extended the test stand by a robot-assisted high-resolution DIC (HR-DIC) system (see Figure <ref>a) <cit.>. It consists of a KUKA iiwa that guides an optical Zeiss 206C stereo light microscope including a Basler a2A5320-23umPro global shutter 16 Megapixel CMOS camera. It enables HR-DIC displacement fields of the crack tip region along the entire fatigue crack propagation. A detailed description of the testing setup including the algorithms for ensuring good image quality and low DIC scattering are presented in <cit.>. Figure <ref>b illustrates the load sequence, including the time slots, during which the HR-DIC data acquisition is conducted.



Before the fatigue crack growth experiment starts, the robot moves the microscope to all positions on the specimen surface, where the fatigue crack is supposed to grow. The area is covered in a checkerboard manner at zero load. At each position a reference image is captured. Automatic algorithms ensure both perfect alignment with the specimen surface and good sharpness and contrast of the images. After this calibration phase, the crack propagation phase starts. The specimen is cyclically loaded until a total crack length of a=27.8mm is measured by a direct current potential drop (DCPD) measuring system (I=100mA, U=60mV). The system stops and the robot moves to the reference image position that best captures the crack tip area. Figure <ref>c  shows the von Mises equivalent strain field. Then 30 load levels are approached, starting with the minimum load up to the maximum load. At each load level, a deformed image is captured. The comparison between the reference image and the deformed images allows to calculate the displacement and strain fields, which is done automatically using the GOM Aramis 2020 software. DIC facet size and distance are 40 × 40 pixels and 30 × 30 pixels, respectively, enabling a spatial resolution of 0.06mm/facet. 



 §.§ Crack opening displacement


To compare the kinematics of crack opening process of the fatigue crack, we use surface displacement data from both the FE simulation and the HR-DIC system. Therefore, the nodal solutions on the surface of the FE model are exported for this study.
We define a COD measurement location, P_cod, within the two-dimensional displacement field by  surface coordinates, P_cod = (x_cod, y_cod), where x_cod is the horizontal distance from the measurement point to the crack tip and y_cod is the vertical distance from the crack path.
Then a second measurement point, P_cod^* = (x_cod, -y_cod), is positioned symmetrically on the opposite side of the crack path, see Figure <ref>a. At these points the displacement component perpendicular to the crack, u_y, is obtained. If a measurement point does not directly coincide with a node or facet, the u_y value is determined by linear interpolation from adjacent points. The local COD value is then calculated as follows:

    COD = 1/2·( u_y(x_cod, y_cod) - u_y(x_cod, -y_cod) )


We plot COD against the load for each discretisation step. The resulting curve, known as the crack opening curve, is shown schematically in Figure <ref>b.



The algorithm for the determination of the crack opening stress intensity factor, K_op, from displacement data follows the methodology recommended by the ASTM E647-15 standard <cit.>. The high quality of the HR-DIC displacement field data allows for a very precise analysis of crack opening. When the measurement points P_cod and P_cod^* are located near the crack tip, three distinct regions within the crack opening curve can be identified (see Figure <ref>b):


	
  (A) A non-linear relationship between load and COD is observed in the lower region, which is characteristic of crack face contact and therefore crack closure. Based on previous research <cit.>, we assume that crack closure will only occur for F < 0.6 · F_max.
	
	
  (B) The middle region shows a linear relationship. The load or stress intensity factor at the crack tip is proportional to the COD, meaning that K ∼COD. This relationship allows the principles of linear-elastic fracture mechanics to be applied.
	
	
  (C) The upper region shows a non-linear relationship caused by large plastic strain accumulation within the plastic zone as reported by Vasco-Olmo <cit.>.	Based on previous research <cit.>, we assume that most plastic strain occurs when loads exceed 80% of the maximum load (F> 0.8 · F_max). 


To determine  K_op, we use the data points from region (B) to fit a line using linear regression. The maximum deviation of a single data point in region (B) from the linear regression model is given by e_max. That means, e_max is a parameter describing the inherent DIC scatter. Then, we define a boundary line by shifting the linear model by 1.2 · e_max to the right, as shown in Figure <ref>b in order to avoid the potential for miscalculations due to the inherent HR-DIC scatter. The first data point that intersects this boundary line, starting from F_min, is considered to be K_op.  
In this paper, we distinguish between three different types of K_op, further explained in Table <ref>. K_pl denotes the first point that leaves the linear-elastic region and crosses the boundary line again.





§ RESULTS AND DISCUSSION
 


We compare the numerical FE results with the  experimental  results obtained from HR-DIC measurements. If there is a satisfactory agreement on the surface, it can be assumed that the crack closure observed in the experiment will behave in a manner similar to that observed in the simulation. This enables the analysis of the 3D aspects of PICC.



 §.§ COD location dependency of K_op,cod

It is widely acknowledged in the fracture mechanics community that the K_op,cod results in DIC displacement field data are sensitive to the measurement location <cit.>. We investigate this dependence through a parameter study of the effects of x_cod and y_cod.  The numerical results are summarized in Figure <ref>, the experimental ones in Figure <ref>. 





The study in Figure <ref> is based on numerical surface displacement fields derived from FE analysis with a crack length of a_f = 27.8mm. Subsequently, the crack tip is fixed to x_ct = 27.8  mm, y_ct = 0  mm.
In order to examine the effect of y_cod, we set x_cod = 0.3mm to be constant. Then, we vary y_cod in a range of y_cod∈ [0.0, 0.6] (Δ y_cod=0.1 mm). Figure <ref>a shows the location of the measurement points P_cod and P_cod^* within the u_y displacement field. The individual measurement points have a distance of Δ y_cod = 0.1mm. The crack opening curves in Figure <ref>b illustrate that the non-linearity due to contact of the crack faces decreases with increasing distance y_cod. In addition, it is observed that the minimum or initial crack opening displacement COD_min increases with increasing y_cod. In contrast to the DIC-based studies by Sutton <cit.>, we observe no dependence on the determined K_op,cod. For the examined position x_cod = 0.3mm we determine K_op,cod = 4.62MPa√(m) = const.
The relationship between the crack opening curves with varying x_cod (Δ x_cod=0.2 mm) is shown in Figure <ref>c and d. In this case, K_op,cod depends on x_cod. Points P_cod close to the crack tip tend to have higher K_op,cod values than those further away. In addition, the evolution of the K_op,cod values shows a gradual, ascending curvature towards the crack tip. Furthermore, the K_op,cod values shown in the crack opening curves  can be divided into two different categories:  Because of R=0.1 and, thus, F_min≠ 0, the crack is only partially closed behind the crack tip due to PICC, schematically illustrated under Figure <ref>d. Furthermore, it seems that especially the section of the closed crack, in which the largest contact pressure between the crack faces exists (see Figure <ref>a), leads to a delayed crack opening near the crack tip. We denote the length of this section behind the crack tip as l_cc. In the numerical study in Figure <ref>, l_cc is known to be 1.36 mm. Figure <ref>d reveals that the K_op,cod form two different cluster depending on if they were measured at P_cod with x_cod<l_cc or x_cod>l_cc. In particular, 
those encircled in black represent K_op,cod at P_cods with direct opposing crack surface contact at F_min (x_cod<l_cc), whereas the dotted highlighted ones indicate K_op,cod at initially open crack surface positions or low contact pressure between both crack faces (x_cod>l_cc). This finding leads to two conclusions: Firstly, with the use of COD measurement it is possible to identify l_cc, if PICC is present, since the K_op,cod based on P_cod with x_cod<l_cc behave differently than with x_cod>l_cc. Secondly, there is still a K_op,cod measurable although the value refers to a P_cod with x_cod>l_cc . This implies the contact close behind the crack tip (<l_cc) still influences the opening behavior of crack sections although their crack faces are not in contact at F=F_min or their contract pressures are so low that they rarely contribute to the delayed crack opening behavior due to PICC. 
We remark that the values of COD_min are occasionally negative. This phenomenon is caused by the FE contact definition and its allowed contact penetration.

Figure <ref> shows the experimental COD measurement results from HR-DIC displacement field data at a crack length of a = 27.8mm. The crack path and the position of the crack tip are identified using the neural networks UNetPath and ParallelNets <cit.> implemented in CrackPy <cit.>, which determined a crack tip position of x_ct = 27.8 mm, y_ct = 0.36 mm. 
In Figure <ref>b we examine the effect of the vertical distance from the crack path y_cod in HR-DIC data. The results reveal no significant deviations in K_op,cod values with variations in y_cod. x_cod=0.6 mm is kept to be constant. This observation confirms that K_op,cod is independent of y_cod, which is consistent with the FE results shown in Figure <ref>b. In addition, the non-linearity of crack closure decreases with the distance from the crack path.
The COD measurement points are located based on the detected crack path at a vertical distance of approximately y_cod≈0.2mm. The horizontal distance between each COD measurement point is Δ x_cod=0.4 mm. 
The results of the x_cod analysis in Figure <ref>d are very similar to the FE results. In particular, the K_op,cod values at points close to the crack tip (x_cod < 1.4mm) are comparable to the FE results, showing the ascending characteristic. In the following, the analysis given indicates that at minimum load F=F_min the crack is closed to l_cc=1.4 mm. For x_cod > 1.4mm the K_op,cod values are almost constant.

In conclusion, our study highlights the negligible influence of y_cod on K_op,cod. Deviations from the results of previous studies <cit.> can likely be attributed to the applied measurement methods. Increasing the distance y_cod results in a flattening of the crack opening curves, and thus increases the effect of measurement noise in HR-DIC data on the determination of K_op,cod. For an accurate characterization of crack closure behavior using HR-DIC, we advise to minimize the noise in DIC measurements and ensure high temporal resolution of the crack opening process. Consequently, despite the negligible influence of y_cod, we recommend to aim for the smallest possible vertical distance of the measurement point P_cod to the crack path. This allows a robust segmentation of the non-linear from the proportional part within the crack opening curve and thus a stable identification of the K_op,cod value. However, both results, numerical and experimental ones, confirm the spatial dependency of K_op,cod along the crack path <cit.> and supports the hypothesis of an increasing K_op,cod with decreasing distance to the crack tip x_cod <cit.>.




 §.§ Identifying plasticity-induced crack closure

As crack surface roughness is not considered, our FE crack propagation simulation only covers the effect of plasticity induced crack closure. Based on the agreement of the K_op,cod behavior in our numerical and experimental studies in Figure <ref>d and <ref>d, we conclude that plasticity induced crack closure is the dominant crack closure mechanism here (AA2024-T3; Δ K = 13.4MPa√(m)). 

In the following, we  aim to  define a criterion to identify PICC that bases on the K_op,cod dependence along the crack path. First, we only take into account for the K_op,cod that bases on P_cod with x_cod < l_cc. Afterwards, we correct the crack opening curves around their initial COD value (COD - COD_min). This is permissible, since we know on the one hand, that contact penetration as allowed within the FE simulation but does not occur in reality. On the other hand, we showed that y_cod, the second reason for a possible horizontal shift of the crack opening curve, does not affect the resulting K_op,cod. Figure <ref>a and c show the crack opening curves of FE and HR-DIC data shifted to a common minimum. Figures <ref>b and d present the  K_op,cod values separately. 



When the COD values of the crack opening curves are aligned to a common minimum, the K_op,cod values exhibit a linear characteristic when related to P_cod with x_cod < l_cc (highlighted in grey in Figures <ref>b and d). Physically, this represents a uniform crack opening towards the crack tip, which is assumed to be a typical characteristic of the crack opening kinematics in presence of PICC <cit.>. Thus, we define the following criterion:

Criterion for plasticity-induced crack closure

Let y_cod be fixed and x_cod,i, i=1,…,n, be x-values of the COD measurement points along the crack path. Let K_op,i be the crack opening stress intensity factor determined for the COD curve at measurement point (x_cod,i, y_cod), and λ_1,...,λ_n be the minimum-shifted crack opening displacements λ_i := COD(x_cod,i,y_cod) - COD_min.
PICC is present and the pre-dominant crack closure mechanism at measurement points I ⊂{1,…,n}, if the following condition holds true: There exists a linear fit, represented by slope and intercept c,k ∈ℝ, such that for all i ∈ I:

    |λ_i - c · K_op,i - k|<ε_picc,

Here, ε_picc > 0 defines the permissible deviation of the points from the linear approximation.

Loosely speaking, PICC is present where K_op,cod∝COD-COD_min with an allowable deviation of ε_picc. This threshold should depend on the discretisation of the crack opening process to reflect the idea that there is less noise in the determination of K_op,cod if the discretization of COD curves is finer, we define ε_picc = Δ K/N_L· s. Here, N_L denotes the number of load steps or sub steps, respectively, used for the discretisation of the crack opening process and Δ K = K_max-K_min. The factor s > 0 is adjustable and can be used to potentially make the criterion more robust. However, it is crucial that the linear character of the K_op,cod values, as illustrated in Figure <ref>, is given. We recommend to keep s as small as possible.
Assuming that the crack continues to open uniformly up to the crack tip, the crack tip opening load K_op,ctod can be estimated by setting the slope of the linear model to zero (c=0) in Equation <ref>. It follows:

    K_op,ctod≈ k


Table <ref> compares the determined values for applying the criterion based on FE analysis with those obtained from HR-DIC data.

Comparing the K_op,ctod result with the sub-step in the FE crack growth simulation when there is initially no contact during crack opening across the crack faces, K_op,FE = 0.43 · K_max = 6.47  MPa√(m), results in a deviation of less than 3 %. The K_op,ctod result based on HR-DIC data is larger compared to the FE data, which can be attributed to significantly higher contact pressures in real fatigue cracks. The high incremental crack advancement in the FE crack growth simulation reduces plastic strain accumulation within the plastic zone, subsequently leading to lower contact pressures on the crack surface. In actual DIC data, measuring crack opening this close behind the crack tip is not feasible due to generally insufficient local resolution of the measurement. However, this method still allows for the determination of K_op,ctod. Generally, the experimental determination of this value is challenging, and a reliable procedure has not been found so far <cit.>. This methodology provides a potential approach to address this issue and thus improve the understanding of crack opening kinematics of plasticity induced crack closure in experimental data.



 §.§ 3D aspects of plasticity-induced crack closure

The FE and HR-DIC results agree on crack opening on the free surface, so we assume that PICC in the experiment behaves similarly to the FE simulation. Contact elements in the crack plane allow us examine the crack face contact through contact pressure distributions. Figure <ref> shows the contact pressure  p_c on the crack surface as contour plots.


All four images indicate that the contact pressure is directed towards the specimen's surface as a result of PICC, aligning with the assertions by Dugdale <cit.> and Newmann <cit.>. They explain this to the significant influence of the plane stress state at the free surface of the specimen. Further FE-based studies have also identified that PICC is particularly prevalent under plane stress conditions <cit.>. Moreover, Figure <ref> shows that, due to the positive load ratio R=0.1, only a limited area directly behind the crack front stays in direct contact. Figure <ref> clarifies that the contact pressure p_c increases with an increasing K_max. The proportion of the contact area A_c to the total area A_ges also increases. At a=57.8mm the entire crack surface behind the crack front is in contact. Additionally, it is observed that the location of the maximum contact pressure p_c,max at a_f=47.8mm shifts from the specimen's surface towards its center. This is particularly evident in Figure <ref>d at a crack length of a_f=57.8mm. 
The variations in the contact pressure distributions can be attributed to the opposing effects of the plastic zone at the specimen's surface and the plane-strain in the specimen's center. As K_max increases, the plane stress state becomes more dominant, resulting in a higher proportion of the crack surface being in contact. Given that contact is primarily focused on the free surface, it is necessary to consider how the external contact of the crack faces influences the crack opening behaviour across the entire specimen thickness. This question is partially answered by Figure <ref>. 
In order to characterise crack opening, the crack opening curves of both CMOD and CTOD nodes (the first row of nodes behind the crack front) along the specimen thickness are examined in Figure <ref>. The basis for this investigation is the FE crack propagation simulation at a=27.8mm from Figure <ref>.
Based on Figure <ref>a, it is apparent that CMOD measurements, regardless of the measurement position in the thickness direction, lead to the same K_op,cmod. Nevertheless, the algorithm-based evaluation (see Section <ref>) detects a K_op,cmod = 2.53MPa√(m), corresponding to 17% of K_max. The surface-near crack closure contact (z=1mm) just behind the crack tip influences the crack opening behavior along the entire specimen thickness. The K_op,ctod values close to the specimen surface feature the largest crack opening values at 47% of K_max. Figure <ref>b, besides illustrating the non-linearity due to crack closure, also shows another non-linear section starting from 70%-80% of the maximum load due to the large plastic strain accumulation in the plastic zone. 

This region is separated by K_pl as introduced in Section <ref> and is determined comparably to K_op. Vasco-Olmo <cit.> identify a similar non-linear behavior in their crack opening curves. They correlate it to the plastic strain accumulation in the plastic zone during crack opening. We support their hypothesis and extend it to three dimensions. The intermediate region between K_op,ctod and K_pl is characterized by linear-elastic crack opening without crack face contact.





 §.§ Effect of crack closure on plastic zone

In the following, we investigate the effect of the crack closure contact on the shape and damage evolution within the plastic zone. We consider as plastic zone (PZ) allplastic  deformations in front of the crack tip. Therefore, we select all elements with plastic deformation after the last opening load step in our FE simulations. The results are illustrated in Figure <ref>. Because of the symmetric FE model, only 1/4 of the entire PZ is displayed. Our results indicate that the PZ exhibits a shape that differs from the dog-bone model <cit.>. Similarly, the numerical shape investigations made by Camas <cit.> and Besel <cit.> are also in agreement with ours, although the effect of crack closure was not considered. Due to the extensive crack advancement, Δ a=0.08mm, employed in the FE simulation, it is only possible to reveal the primary PZ. Consequently, the cyclic PZ cannot be evaluated.
In order to analyse the impact of crack surface contact, we conducted FE crack propagation simulations with and without contact definitions. Subsequently, we compared the shapes of the resulting primary PZ. The use of elements with a cuboid shape facilitates quantitative assessments concerning the shape characteristics of the primary PZ. Figure <ref> shows the PZ shapes with active contact definition at four different crack lengths. Table <ref> lists the results of selected descriptors for the shape of PZ. A summary about the parameter is illustrated in Figure <ref>a.



Comparing the descriptors listed in the Table <ref>, the deviations are all less than the size of one element, i.e. less than 2 % in all analysed parameters. Therefore, we conclude that the influence of the crack closure contact on the shape of the primary PZ is negligible. Whether the contact is confined to the crack surface near the free surface or encompasses the entire crack surface appears to be unimportant. Furthermore, the contact pressure distributions correlate with the shape of the primary plastic zone, as shown in Figure <ref>. Figure <ref>b plots the fraction of crack surface being in contact s_c to the entire specimen thickness (see Figure <ref>a). The evaluation is based on the second row of elements behind the crack front to avoid numerical singularity effects that can result from the crack propagation algorithm. As illustrated in Figure <ref>b, full-crack surface contact is present when the plastic zone is greater in width than in height. Since the shape of the plastic zone depends on both K_max and the specimen thickness t <cit.>, this statement is so far only valid for the specimen under investigation (t=2mm, AA2024-T3). Nevertheless, the impact of crack closure contact on the shape of the cyclic plastic zone remains to be investigated.


Figure <ref>b compares the shape of the plastic zone at the surface based on FE results with those based on the HR-DIC analysis at a crack length of a=27.8mm. Both plastic zones are delineated from the linear-elastically deformed material by the yield strength R_p0.2=350MPa. It is observed that both shapes almost coincide indicating that the experimental results agree with the numerical findings. Furthermore, since the upper wing of the plastic zone is greater in height than in width, it indicates that the contact is only present near the free surface of the specimen. We neglect the lower wing because it has been found that its shape is strongly influenced by a HR-DIC measurement artifact resulting from low quality speckle pattern. Nevertheless, the plastic zone based on HR-DIC data is estimated to be slightly larger, attributed to the inherent measurement noise in HR-DIC data.
Because of the good agreement between the numerical and experimental results, we can analyse damage in the plastic zone due to plasticity induced crack closure. Therefore, we use the accumulation of plastic energy over a load cycle as damage describing parameter, following the suggestion by Vormwald <cit.>. In our analysis, we compare the results of an FE simulation with contact definition to one without contact definition. This comparison helps to identify the influence of the contact on the energy development within the plastic zone. Figure <ref> illustrates the distribution of plastic energy d U_pl/d N along the crack front z.



The two curves shown in  Figure <ref>a deviate less than 1%, but d U_pl/d N with contact definition is slightly larger. This effect results from the numerical singularity that arises from the first FE element behind the crack tip and whose influence is slightly increased by the contact pressure. We observe the same trend for the crack lengths shown in Figures <ref> b, c and d. The increased portion of the crack face being in contact combined with larger contact pressures increases the numerical singularity stress within the crack front element. This causes that the distance between the with contact and no contact curves increases. Only for larger crack lengths a_f=57.8mm the contact pressure appears to reduce the accumulation of plastic energy per load cycle. The contact pressure is no longer focused towards the free specimen surface.  At a=57.8mm a plane stress state is pre-dominant for the entire plastic zone and the crack surface behind the crack front is in full contact. The findings indicate that when the contact pressure is concentrated on the free specimen surface, there is no notable impact on the plastic energy accumulation from PICC. However, when the maximum contact pressure point is shifted towards the specimen center and the entire crack face is in contact, PICC appears to exert a reducing effect on the plastic energy accumulation within the PZ.



 §.§ Fracture surface analysis

Finally, we examined the fracture surface by scanning electron microscopy at different scales to find indications of crack surface contact near the surface. Figure <ref>a shows the fracture surface section whose displacement field is shown in Figure <ref>c and to which all COD analyses in Figure <ref> refer. The lower edge of the fracture represents the free specimen surface that has been investigated by the robot-supported HR-DIC measurement system. Before a=27.8mm, two different damaging mechanisms can be identified within the fracture surface. In the center, as shown in Figure <ref>b, striations are observed. This aligns with the literature towards intrinsic fatigue crack growth mechanisms in AA2024 <cit.>.


However, since the previous analysis has shown that the crack closure contact due to PICC is focused towards the crack surface edge, we study the intrinsic mechanisms in that particular region. Figure <ref>c shows that the fatigue crack propagates in that region by ductile shear failure mechanisms leading to a rough fracture surface close to the specimen surface. At Δ K=13.4MPa√(m), fatigue crack growth in t=2mm thin MT(160) specimens made of AA2024-T3 tends to develop shear lips. It follows that the surface roughness in that region is increased. However, since the crack opening kinematics are still identical to  those observed in the FE crack propagation simulation, this indicates that the increased surface roughness due to the shear lips has a minor influence on the pre-dominant PICC mechanism in this section of AA2024-T3 fatigue cracks. That finding supports research performed by Materna <cit.>. They concluded, based on numerical studies, that fracture surface roughness does not affect the crack closure behavior if PICC is the pre-dominant mechanism. 





§ CONCLUSIONS


In this study, we investigate the crack closure behaviour of AA2024-T3 based on numerical simulations and experimental data. We find an excellent agreement between the crack opening curves based on FE and HR-DIC data. Given that the FE crack propagation model only considers the effect of PICC, we conclude that PICC is the dominant crack closure mechanism in AA2024-T3 fatigue cracks in our studies for an MT(160) specimen with t=2mm and a Δ K=13.4MPa√(m). In addition, the agreements between simulation and experiment leads us to the following conclusions:




	
  * PICC can be identified by using the K_op,cod dependence on the COD measurement location. Here, K_op,cod only depends on the horizontal distance to the crack tip position (x_cod) and is independent of the vertical distance. However, we recommend to position the measurement points as close as possible to the crack path enabling a more stable identification of K_op,cod.
	
	
  * Using the variation of K_op,cod values along the crack path, we introduce a new method for determining the opening value directly behind the crack tip K_op,ctod. If PICC is present and the crack opening curves are shifted to a common minimum, the K_op,cod values form a linear relationship, that can be approximated by linear regression.
	
	
  * The study supports the general assumption that PICC induces crack surface contact focused towards the free specimen surface. However, we show that even if there is only crack surface contact near the free surface of the specimen, it influences the crack opening behaviour throughout the specimen thickness.
	
	
  * Considering a fatigue crack in AA2024 and a t=2mm thick MT(160) specimen, we found that the shape of plastic zone  is correlated to the crack surface contact situation. If the width of the plastic zone is larger than its height, this indicates that the entire crack surface close to the crack front is in direct contact. 
	
	
  * Regarding the retardation effect of PICC on the damage accumulation within the plastic zone, we find  that the crack closure contact does not affect the plastic energy accumulation if the crack closure contact is directed towards the free surface of the specimen. If the entire crack surface is in contact, this seems to reduce the energy accumulation. 
	
	
  * Based on fractographic investigations, we observe that fatigue cracks in AA2024-T3 at Δ K=13.4MPa√(m) tend to develop shear lips. However, the increased surface roughness does not appear to affect the PICC mechanism. 





§ ACKNOWLEDGEMENTS



We acknowledge the financial support of the DLR-Directorate Aeronautics. This work was supported by the Deutsche Forschungsgemeinschaft, Germany (DFG) via the project Experimental analysis and phase-field modeling of the interaction between plastic zone and fatigue crack growth in ductile materials under complex loading (grant number BR 6259/2-1). Furthermore, funding came from the Federal Ministry for Economic Affairs and Climate Action, Germany on the basis of a decision by the German Bundestag, within the framework of the aerospace program LuFo-VI of the project "Intelligent FSW Process Monitoring" (Funding ID 20W2201E).




§ DATA AVAILABILITY

The apdl finite element simulation code and data as well as the experimental displacement fields will be publicly available on Github (<https://github.com/dlr-wf>) and Zenodo (<10.5281/zenodo.13643861>).



§ COMPETING INTERESTS

The Authors declare no Competing Financial or Non-Financial Interests.



§ AUTHOR CONTRIBUTIONS

F.P. conceived the idea, conducted the simulations and experiments, evaluated the results, joined the discussions and wrote the manuscript. D.M. joined the discussions and wrote the manuscript. E.B. joined the discussions and wrote the manuscript.



