Three-dimensional simulation of crack propagation in the central-cracked flat and stiffened plates under uniform tension

In this paper, the crack propagation behavior in flat and stiffened plates with central-cracked is studied based on the theory of linear elastic fracture mechanics and 3D finite-element method. The magnitude and distribution of the stress intensity factor in a 3D stiffened plates are affected by the out-of-plane bending and loading modes. Initially, for validating this method, the behavior of crack propagation in a central-cracked flat plate (unstiffened), followed by the propagation behavior and the crack arrest effects on stiffened plates by stiffeners, as well as the out-of-plane bending effect on the geometry correction function distribution are studied; However , the results are compared with the results of the referenced article. In order to analyze the effect of stiffeners in preventing crack expansion, stiffeners with variable thickness and height were used. It should also be noted that the crack propagation behavior in the stiffened plate is analyzed in two methods. In the first analysis, after the crack tip reaches to the stiffeners, they have not fracture and the crack only expands in the stiffened plate (in the stiffened plate type 1 and 2). In the second analysis, with the crack growth in the stiffened plate and the crack tip reaching to the stiffeners, Also they fractured and cracks in the stiffened plate and stiffeners are propagated (only in the stiffened plate type 2). In both analyses, the magnitudes of the geometry correction function in the middle-plane, the crack arrest effect by the stiffeners (Δβ s ), and the out-of-plane bending effect (Δβ) were obtained. Then, the changes in these magnitudes were investigated. It was observed that the stiffeners thickness had a negligible effect on the distribution of the geometry correction function and, in general, the stiffeners had a significant role in preventing the crack growth compared to the flat plate, and the more the crack tip is closer to the adjoining stiffeners; the stiffeners effect in the crack arrest will increases. Furthermore, with attention to the geometry correction function that is studied for different crack sizes, it was determined that the geometry correction function and stress intensity coefficient with increasing the length of the crack in the plate, increase.


Introduction
Stiffened plates are important structure components that generally are manufactured from metal stiffened plates with stiffeners and designed to withstand different loading conditions. These plates are widely used in the industry due to light-weight and high strength and stiffness. However, these structures are often fractured due to the crack propagation in plates or stiffeners. As a result, the crack growth behavior in the stiffened plate is more complicated than the flat plate and has been studied by many researchers [1]. In recent years, linear elastic fracture mechanics have been used extensively to predict fracture in brittle materials, and one of the basic parameters in this case is the stress intensity factor used to determine crack growth as well as the crack critical length. The coefficient of non-dimensional intensity stress, which in some cases being called the beta factor, is called geometry correction function. This factor is used to calculate geometry effects and also usually for fracture analysis. In some references, the effects of length and width on stress intensity factor in plates with crack were analyzed, and geometry correction function was calculated [2]. The stress intensity factors in stiffened plates with crack are influenced by the type of stiffener, its size and location, and these effects should be calculated in the process of crack propagation analysis. In the earlier works, the boundary displacement method was used to derive solutions for stress intensity factors [3]. The comparison of damage tolerance for various stiffened plates was performed by Nesterenko [4] and Zhang et al. [5] with the finite-element method (FEM) and experiments. Won et al. [6], Garcia et al [7] have obtained the characteristics of three-dimensional stress fields in central-cracked flat plates. Chang et al., Using solid elements, investigated the stress intensity factors of skin-stiffener with inclined cracks [8]. In this paper, using the linear elastic fracture mechanics, theoretically the parameters, including the geometry correction function components, the non-dimension geometry correction function and the intensity stress factor are investigated. Then the propagation behavior and crack arrest effect of the central-cracked flat and stiffened plate under uniform tension are studied using the ABAQUS finite-element software. Also, the magnitude of the geometry correction function in the middle-plane, the arrest crack effect by the stiffeners and the out-of-plane bending effect at stiffeners different heights, were obtained with two different analysis methods. In the first analysis, the crack growth was performed only on the stiffened plate, and the stiffeners were not fractured (on the stiffened plate type 1 and 2). In the second analysis, the crack in the stiffened plates and stiffeners has been expanded, and also the stiffeners have been fractured (only in the stiffened plate type 2). Finally, the analysis of the geometry correction function magnitude in the middle-plane, the arrest crack effect by stiffeners and the out-of-plane bending effect are discussed. The results of this work have been compared with the results of Huang et al. [9].

Investigating central-cracked structures
In this paper, three types of central-cracked structures are investigated. First, the flat plate with a through-the-thickness crack was examined by three different thickness of the plate along different cracks. In the following, similar to the flat plate, the study of the stiffened plates with central-cracked has been studied. The flat plate and its coordinate system are shown in Fig. 1. The geometry parameters for the flat plate, including H/W=1 and t/W=0.015, 0.025, 0.055, and normalized crack's length are specified in Table 1. Also, the two ends of the flat plate are under uniform tensile stress equal 120MPa (σ = 120MPa).  The stiffened plate type 1 is including a plate and stiffeners that are symmetrical about the central plane (yoz) of the plate (Fig. 2). The two upper and lower ends of the plate, under uniform tensile stresses, have located that equivalent to 120MPa (σ=120MPa), and the two ends of the stiffeners are assumed to be free. Also, the geometry parameters of the stiffened plate similar to the flat plate, such as H/W=1, and t/W=0.015, were considered, but stiffeners were investigated in different sizes (height hs and thickness ts). In the investigation the stiffened plate, at first two same-size stiffeners are located on the two edges of the plate, namely, the value of 2ds/W (ds is half the stiffener spacing from each other) is approximately equal to 1. In this stiffened plate, the crack is propagated only in the plate, and the stiffeners do not fracture. In stiffened plate type 1, flat plate and stiffened plate are assumed homogeneous and isotropic with the Poisson's ratio ν = 0.3 and the young's modulus E=200GPa, and on the other hand, the young's modulus of the stiffeners is also considered equal to the young's modulus of the plate. In the following, the crack propagation behavior in the stiffened plate type 2 has been studied with uniform tensile stresses. Stiffened plate type 2 consists of a plate and two same-size stiffeners, such that the upper and lower ends of the plate are under tensile stress, and stiffeners are assumed free in two ends. Furthermore, the stiffeners are symmetrical about the central plane (yoz) of the plate and are connected to the plate using the tie constraint. The spacing of the stiffeners to the symmetry plate (ds) is equal to 0.55 half-width of the plate (W/2) and the spacing of the stiffeners from each other (2ds) is 1/1 half-width of the plate (W/2) (Fig. 3). The geometry parameters for the stiffened plate, including H/W=1 (H, W are respectively plate height and plate width) and t/W=0.018 (t is plate thickness). Also, magnitudes of the normalized crack's length are shown in Table 2. The stiffeners thickness (ts) is assumed equal to the plate thickness (t) and the stiffeners are investigated only at different heights (the studied hs/t is 5.5, 10.75, and 13). As noted above, in stiffened plate type 2, at first with crack propagation in the plate, stiffeners was not fractured and the crack propagation occurred only in the plate, and then, with the crack growth and the crack tip reaches the stiffeners, they also fractured and the cracks in stiffened plate and stiffeners have expanded. The stiffened plate, homogeneous and isotropic is assumed to be with the young's modulus E=205GPa and the Poisson's ratio ν=0.3.

The theoretical formulation of cracked structures
In linear elastic fracture mechanics, stress intensity factor is an important parameter for displaying the stress field at the crack tip. For central-cracked flat plate that under the uniform uniaxial tensile stresses ( Fig. 1), non-dimensional stress intensity factor is calculated from the following relation: In this relation, K is the mode-I stress intensity factor, σ is the applied stress, a is half the crack's length of the plate, and β is the geometry correction function. Using the results of the twodimensional correction function for the flat plate, a simple equation for extracting the threedimensional geometry correction function is presented [16]: In this relation, 3 is the geometry correction function at the middle-plane of a Threedimensional flat plate, 2 is the geometry correction function of the two-dimensional flat plate and ν is the Poisson's ratio. It should be noted that the two-dimensional analysis does not usually take into account the stress intensity factor in thin plates.

Fig. 4. Decomposition of a stiffened plate
The linear fracture analysis of the stiffened plate is based on the flat plate analysis. As shown in Fig. 4, for investigation of the stress field, the stiffened plates are divided to the plate and the stiffeners and also the transverse shear force between the plate and the stiffeners is ignored. Using the principle of superposition, the stress intensity factor for stiffened plate is defined as follows.
In this relation, Q is the longitudinal shear force that applied to the plate by stiffeners. KQ is the corresponding stress intensity factor applied by Q and generally ≥ 0, βp is the geometry correction of the flat plate, and a is half the crack's length in the plate. In the upper relation, in the right-hand side of the relation, the stress intensity factor of the flat plate is correspond. When the central-cracked stiffened plate is under uniform tensile stress, the out-of-plane bending occurs as a result of geometry asymmetry (as shown in Fig. 2). Due to out-of-plane bending, a part of the cross-section (xoz) is placed under additional tensile stress and applied to the other side of the compression stress. Additional compression and tensile stress affect the magnitude and distribution of the stress intensity factors, thus the crack propagation effect on the stiffened plate will be different. So it's best to be investigated the effect of out-of-plane bending. With regard to the effect of out-of-plate bending, relation 3 can be written as follows: In the above relations, M is the additional torque caused by the out-of-plane bending, σM is the additional stress caused by the out-of-plane bending, Ix is the moment of inertia about the X-axis for the stiffened plate, and ́ the distance from the calculated location to the neutral plane of the stiffened plate. In order to simplify the problem, and considering that the stress intensity factor changes with the stiffened plate thickness direction, the middle-plane (z/t = 0.5) of the stiffened plate was considered as a neutral plane. In this case, the stress intensity factor on the middleplane (Kmid) and the maximum stress intensity factor (Kmax) is: Where σMmax is the maximum additional tensile stress caused by out-of-plane bending, and zmax is the distance of the neutral plane to the cross-section where the additional tensile stress reaches its maximum. The stress intensity factor can be converted to non-dimension parameters, so the geometry correction function in the middle-plane (βmid) and the maximum geometry correction function (βmax) in the stiffened plate is equal to: In this paper, the discussion and analysis of the crack arrest effect and its propagation are limited to the geometry correction function. In order to determine the role of stiffeners in 3D conditions, the crack arrest effect on stiffeners and 3D stiffened plate is proposed as follows [16]. The crack arrest effect by stiffeners (Δβs), which is the same influence of -Q on the geometry correction function in the plate, is equal to: The parameters of this relation are already defined. To illustrate the effect of the torque caused by out-of-plane bending on the stress intensity factor, the following parameter is defined: Where Δβ is the increase of the geometry correction function from the middle-plane to the maximum value. Compared to the flat plate, the crack arrest effect in the stiffened plate (Δβsp) can be defined as: Relation (13) parameters have already been defined. The concept of the crack arrest effect is based on the 3D geometry correction function. If the values of βp, βmax, and βmid are extracted from finite-element analysis, the values of Δβ, Δβs, Δβsp can be obtained.

Simulation and verification of the finite-element method
The stiffened plate model in ABAQUS finite-element software version 6.14.3 has been used and in the simulation of crack propagation behavior, Contour Integral method has been used. To calculate the stress intensity factors along the crack tip in a three-dimensional stress field, the domain integral function proposed by Shih et al. [10] has been used and the mode-I stress intensity factors are calculated from the following relation: In this relation, J is the domain integral; E is the young's modulus, and ν is the Poisson's ratio. In the article, three-dimensional brick elements are used, which are suitable for examining the stress intensity factor along the plate length direction. In this paper, due to the symmetry of the flat plate, only half of it was modeled with 8-node linear brick and reduced integration element (R8D3C) in ABAQUS software. Furthermore, three types of meshes, with 185500 elements, 275000 elements and 555000 elements all of which are of the type (R8D3C), were used to examine the convergence of finite-element analysis. In the end, the difference in the values of geometry correction functions on the middle-plane was estimated to be less than %0.106. It should be noted that the geometry correction functions obtained in Table 1 are calculated with regard to 185500 elements of the mesh. Similarly, half of the stiffened plates (type 1 and 2) were modeled, and a mesh with a number of 228000 elements (R8D3C) was used in the finite-element analysis. The number of mentioned components varied slightly depending on the size of the stiffeners and the cracks. Three samples of the finite-element models are shown in Figures (5) and (6). It should also be noted that in the stiffened plate analysis of the type 2, considering that the strain field is focused on the crack tip, and in the present work, the analysis of the middleplane (z/t = 0.5), in order to increase the accuracy Integral Counter in the calculation of stress and strain, singularity was used. In order to validate the finite-element method, the crack propagation effect and crack arrest effect on the flat plate for crack different sizes and in three thickness of the flat plate is investigated. In the table (3), using the values of the domain integral function and the stress intensity factor, the geometry correction function in the middle-plane (βmid) was calculated with the finite-element method (Table 1). Then its results were compared with the referenced-article results and subsequent the error percentage was calculated. In the following, the geometry correction function figure along plate thickness direction for crack different sizes and in three thickness of the flat plate is investigated. And, for example, the figure has derived for t/W = 0.025 (Fig. 7). Furthermore, the values of the geometry correction function were plotted in different crack sizes for t= 0.015, 0.025, 0.055 (Fig. 8) and compared with the figures of the referenced-article.

Analysis of geometry correction functions related to flat plates
The geometry correction functions are not uniform along plate thickness direction and reach at z/t=0.5 to the maximum value. Furthermore, as shown in Fig. 7, the distribution of geometry correction functions on both sides of the middle-plane is symmetric. Therefore, it can be concluded that the process of changing the geometry correction function along the crack is different and its distribution symmetrical rather to the middle-plane. This symmetry is observed in relation to the middle-plane, which has the maximum value of the geometry correction function, in different thicknesses of the flat plate as well as in different crack sizes. By observing Fig. 8, it can be said that the changes in the geometry correction function relative to the different sizes of the crack in the different thickness of the plate follow a similar process, such that, as long as the crack length increases, the value of the geometry correction function also increases. Furthermore, the process of changing the geometry correction function in the referenced-article is presented in Fig. 8 and as it is shown, the process of changes in the geometry function in different sizes of crack in the present work and the referenced-article is the same and the values of the calculated geometry correction function are also acceptable percentages in comparison to the referenced-article.
The values of the geometry correction function at the middle-plane in the three-dimensional analysis are greater than the two-dimensional analysis, which indicates the low accuracy of the two-dimensional analysis in thin plates. According to the previous study [12], and comparing the results of the finite-element method with the results of relation (2), it is determined that with increasing plate thickness, the results of the three-dimensional analysis are close to twodimensional. These results are consistent with the referenced-article and confirm the validity of the three-dimensional finite-element method for crack analysis.

The crack propagation behavior in the central-cracked stiffened plate has two symmetric stiffeners at the edge (stiffened plate type 1)
Because of the asymmetry of the stiffened plate compared to the middle-plane (z/t = 0.5), when the plate is under uniform tensile, the out-of-plane bending occurs (Fig. 9). Therefore, as shown in Fig. 10, the distribution of the geometry correction function is not symmetrical relative to the middle-plane. It is worth noting that the out-of-plane bending has changed by the change in the stiffened size (height hs and thickness ts), which affects the distribution of the geometry correction function.
For a more precise examination of the out-of-plane bending effect, the distribution of the geometry correction function was investigated at different heights and stiffeners thickness. Initially, the stiffeners thickness value was fixed assumed, and the distribution of the geometry h s /t =2 t s /t = 1 correction function along plate thickness direction at stiffeners different heights was examined (Fig. 11). In Fig. 12, keeping stiffener height )hs) as constant, the distributing processes of the geometry correction function in stiffeners different thickness were investigated. It is observed that the slope variations of the geometry correction function distribution curves are low and, as a result, the stiffener thickness has a slight effect on the geometry correction function distribution. In Figures 10 to 12, there are two important points in the geometry correction function distribution curve, namely, βmax, and βmid, which can be evaluated by the out-of-plane bending effect with respect to these values in the figures of Δβ (Figures 13 and 14). With respect to the changes in the Δβ Curves in Figures 13 and 14, it can be seen that Δβ has changed from 0 to %10 with a change in stiffeners height and from 0 to %1.2 with a change in stiffeners thickness. This indicates the stiffeners height effect on the geometry correction function distribution and the limiting effect of the stiffener thickness on the function distribution. In fact, with the increasing effect of the stiffener height on the Inertia moment about the X-axis (Ix) in comparison with the stiffened thickness, this fact can be seen. The values of βmid in different heights (Fig. 15) and various stiffeners thicknesses (Fig. 16) are shown. Also, the corresponding Δ s, namely the influences of -Q on the geometry correction function, are indicated in Fig. 17 and 18. From Figures 15 through 17 it can be concluded that the crack arrest effect by the stiffeners depends on the distance the crack tip from the adjoining stiffeners. When the normalization length is less than 2a/W=0.45, there is slight difference between the flat plate and the stiffened plate, but when the normalization length is increased from 2a/W=0.55, the stiffeners prevent crack propagation. In other words, the crack arrest effect by the stiffeners increasing with its size increases and normalized crack's length more than 2a/W=0.55 also intensifies. Because of the stiffeners are located at the plate edges, and only the crack is propagated in the plate and compared to the flat plate, it can generally be said that the stiffeners effectively prevents the crack's growth. This problem is along normalized crack more than 2a/W=0.55, and when the crack tip approaches to its adjoining stiffeners, can be seen.

Analyses of the crack propagation behavior in central-cracked stiffened plate under uniform tensile (stiffened plate type 2)
In stiffened plate type 2, as long as the crack tip does not reach the stiffeners, it can be said that the crack arrest effect by stiffeners is similar to the stiffeners plate type 1 (two stiffeners on the edges of the plate). With the crack propagation on the stiffened plate type 2 and the crack tip reaches to the stiffeners; the crack propagation behavior has changed and takes a different shape. Regarding Figures 20 and 21, it can be seen that changes in the size and location of the stiffeners are affected by the crack arrest. As shown in Fig. 20, with the increase of hs/t values, the magnitudes of the geometry correction function in the middle-plane decreased. This decrease is observed in the βmid value by increasing the hs/t value in both examined methods (Fig. 19). Also, as shown in Fig. 20, before the crack tip reaches to the stiffeners, the βmid values (crack growth only in the plate) have been less than the βmid values in the stiffened plate with two equal stiffeners placed on the plate edge. This indicates the more increase in the crack arrest effect by stiffeners with approaching them from the plate edges towards the symmetry-plane. Furthermore, the results were compared to the results of Huang et al., and good conformity was observed. In another analysis, with the crack growth on the plate and the crack tip reaches to the stiffeners, they also are broken (parts 3 and 4, Figure 19). According to Fig. 20, when the stiffeners have been broken, the βmid values are only slightly reduced and increased by the crack growth on the plate and stiffeners. As the cracks continue to propagation, the stiffeners are completely broken, consequence, the crack arrest effect is destroyed by the stiffeners and the βmid behavior approaches to the flat plate. According to the two analyzes mentioned above, it can be seen that the arrest crack effect by stiffeners was much higher when the crack growth was performed only in the plate, and the stiffeners were not broken. The out-of-plane bending effect can also be analyzed using Fig. 22 and examining the Δβ values. With the crack growth in the plate and before the crack tip reaches to the stiffeners, as can be seen, the out-of-plane bending effect is negligible. By continuing the crack growth in the plate and the crack tip reaches to the stiffeners, the out-of-plane bending effect in the two analyzes noted (Fig. 19) has increased sharply. This increase indicates the importance of the out-of-plane bending effect, especially after the crack tip reaches to the stiffeners, and rejects the assumption of ignoring the out-of-plane bending effect. The results obtained in the two analyze mentioned above are compared with the results of Huang et al., And the proper conformity has confirmed the validity of the results [14].

Conclusion
In this paper, the crack propagation behavior in flat and stiffened plates with central-cracked is studied based on the theory of linear elastic fracture mechanics and 3D finite-element method.
The magnitude and distribution of the stress intensity factor in a 3D stiffened plates are affected by the out-of-plane bending and loading modes. The obtained results are as below: 1. Stiffeners have a significant role in preventing crack growth in stiffened plates compared to flat plates. 2. With the crack tip reaches to the adjacent stiffeners; the stiffeners effect will increase in the arrest crack.
3. According to the study of stiffened plates at different heights and thickness of the stiffeners, it was observed that the stiffeners thickness had a negligible effect on the geometry correction function distribution. 4. The geometry correction function and the stress intensity factor increases with increasing the crack length in the plate. 5. The out-of-plane bending effect before the crack tip reaches to the stiffeners is negligible and its importance, especially after the crack tip reaches to the stiffeners, has increased (it has increased significantly), and rejects the assumption of ignoring its. 6. In stiffened plate type 2, if the stiffeners are completely broken, the crack arrest effect by the stiffeners is eliminated and the βmid behavior approaches to the flat plate.