Damage evolution in plates subjected to fatigue loading

The aim of this paper is the numerical prediction of the cracking path followed by a surface crack front in plates constituted of different materials (determined by the exponent m of the Paris law), subjected to cyclic tension or cyclic bending loading. To this end, a numerical modelling was developed on the basis of the discretization of the crack front (characterized with elliptical shape) and the crack advance at each point perpendicular to such a front, according to a Paris law, using the stress intensity factor (SIF) calculated by Newman and Raju. Results show that the crack leads to a preferential propagation path that corresponds to a very shallow initial crack with a quasi-circular crack front. The increase of the Paris exponent produces a quicker convergence during fatigue crack propagation from the different initial crack shapes.


Introduction
The plate, a simple geometry of structural application, often presents cracks on its surface that can easily propagate due to fatigue, a reason why studying them is of great interest in fracture mechanics. The stress intensity factor (SIF) in a superficially cracked plate, subjected to tension or bending loads, has been calculated by several methods [1][2][3][4][5][6][7][8][9][10][11]. For tension loading the SIF increases with crack depth, showing its maximum value at the centre of the front for semi-elliptical cracks (with a ratio between semiaxes of 0.2) and at the point of the surface for semicircular cracks [1,2]. For bending moment the SIF is highest at the crack surface from relative depths greater than or equal to 0.6, regardless of the geometry of the crack front [1].

Numerical modelling
A computer application in Java programming language was developed to study the propagation path of a symmetric surface crack on the cross section of a plate (figure 1), subjected to fatigue (cyclic) loading in the form of remote uniform tension or pure bending moment (under Mode-I loading conditions).
The SIFs used were those obtained by Newman and Raju [1][2][3] by a 3D finite element analysis and the nodal-force method. They adjusted their results to an equation [3] that calculates the dimensionless SIF with parameters a/b, a/t, b/w and  characterizing each point of the crack front (figure 2). Although other equations do exist, this one has been the most widely used in the bibliography due to its excellent results. This equation (valid for 0 ≤  ≤ π, b/w < 0.5, 0 ≤ a/b ≤ 2.0 and a/t < 1.25(a/b + 0.6) if 0 ≤ a/b ≤ 0.2 or a/t < 1 if a/b > 0.20, but for a/t > 0.8 the accuracy has not been well established) has been tested experimentally and numerically by various methods [1,5,8,9,[12][13][14]16,20]: measurements of the SIF K (photoelastic and numerical), the geometrical evolution of the front, the crack growth rate, the failure stress, etc. The SIF expression has the following form: wherein S t is the remote tensile stress, and S b the stress caused by the moment (maximum remote bending stress), The parameter F s is obtained by the following equation, The parameters needed to obtain F s and H s (M 1 , M 2 , M 3 , g, f  , f w , H 1 , H 2 , p and ), as well as the parameter Q, are calculated by means of the mathematical expressions given in the classical paper by Newman and Raju [3]. Such parameters are, in general terms, functions of the geometry of the cracked plate, of the kind of loading (tensile remote stress or bending moment) and of the specific point at the crack front P (cf. figure 2). The basic hypothesis of this modelling is the assumption that the crack front can be characterized as a semi-ellipse with its centre on the plate's surface and that fatigue propagation takes place in a direction perpendicular to such a front, following the Paris-Erdogan law [19], where C and m are the Paris parameters.
To discretize the crack front, modelled in a semi-elliptical way, the front was divided into z segments of equal length, using Simpson's rule. The point on the surface of the plate, where a plane stress state does exist, was not taken into account. Subsequently, each of the points (i) was moved perpendicularly to the crack front, in accordance with the Paris-Erdogan law, so that the maximum crack increase, Δa(max), associated with the point of the maximum SIF range, ΔK(max), was kept constant.
From this maximum increase, and using the ratio of the SIF range at a given point of the crack to the maximum SIF range over the crack front, the advance of each one of the front points (Δa i ) can be obtained as follows, Keeping in mind Newman and Raju's equations [3], the advancement at each point of the crack front can be calculated using the expression, considering that H s is equal to 1 in the case of tension loading.
The new points, fitted by the method of least squares, form a new semi-ellipse, so that the crack advancement progress is repeated iteratively, according to the flowchart shown in figure 3. Due to the existing symmetry, only half of the problem was used in the calculations.

Numerical results
A convergence study was carried out in order to determine the number of parts in which each semiellipse is divided (crack front) and the maximum crack depth in the iterations. The parameters selected for the calculations were z = 24 and Δa(max) = 0.00001t. Figure 4 plots the evolution of the aspect ratio a/b with relative crack depth a/t. It shows how cracks propagate from distinctive initial geometry -relative crack depth (a/t) 0 = {0.02, 0.1, 0.2, 0.3, 0.4, 0.5} and aspect ratio (a/b) 0 = {0.2, 0.5, 1, 1.5, 2}-in plates of different materials (m = 2, 3 and 4), subjected to tension loading or bending moment.
The crack advance from different initial geometries of the surface defect (a/b-a/t curves) tends towards a preferential propagation path, which corresponds to that of an initially very shallow crack with a quasi-circular front (a/b ~ 1).