Numerical Study on Ductile Failure Behaviours of Steel Structures under Quasi-Static Punch Loading

<jats:p>A reliable finite element procedure to simulate shear-dominated ductile fractures in large-scale, thin-walled steel structures is still evolving primarily due to the challenges in determining the failure criterion of metal materials under complex stress states. This paper aims to examine the accuracy of the modified Gurson–Tvergaard–Needleman (GTN) model considering the shear failure in simulating the ductile fracture of steel plate structures under quasi-static punch loading. The modified GTN damage models are performed by the ABAQUS user-defined material subroutine (VUMAT). The void-related parameters and shear damage parameters of Xue’s and the N-H modified GTN models are calibrated from test specimens with various geometries corresponding to different stress triaxiality and shearing conditions. The damage evolution associated with shearing of voids in the modified GTN models has strong influences on the stress triaxiality versus plastic strain under complex stress states, especially for the shear-dominated loading conditions. Based on the original GTN model, Xue’s and the N-H modified GTN model with calibrated material parameters, a numerical comparative study examines the ductile fracture of steel non-stiffened plates and stiffened plates under punch loading. Benchmarked against the experimental studies, the numerical simulations demonstrate that the shear-driven void evolution in the modified GTN model imposes significant effects on the load–displacement responses as well as the onset and extension of ductile fractures in steel plates under punch actions. The N-H modified model with calibrated shear damage parameters shows a better correlation with the ductile fractures in steel plates observed in the experiment than the original GTN model and Xue’s modified GTN model. As a result of this study, the modified GTN model considering shear action can be applied for practical applications in the crashworthiness assessment of ship collision and grounding.</jats:p>

section_content: Ship collision and grounding accidents have become common events with increasing tonnage and busier traffic routes. In a collision or grounding accident, the penetration and damage of the ship structure may cause loss of life and property as well as potential pollution [1]. In most cases, ship collision or grounding events are the scenarios of structures subjected to local quasi-static loads with a low loading speed and large loading mass. The crashworthiness analysis of a structure subjected to localised quasi-static punch load is widely used in the safety assessment in the early design stage of a ship. It is very important to develop an effective and accurate assessment procedure for the crashworthiness of a structure under localised load. Consequently, a comprehensive understanding of the failure mechanisms, including the fracture failure driven by shear actions, in marine structures under local quasi-static loads thus becomes essential to enhance the integrity of these structures. 

The explicit finite element method is often used to analyse accidental responses of large-scale and thin-walled steel structures due to the rapid development of computational technology. However, the simulation of cracking evolution in marine structures is a challenging task in the ship collision or grounding analysis, especially in the determination of the fracture criterion under large shear actions. Numerical simulations for the penetration responses of a structure including highly nonlinear elastic-plastic deformations until rupture are also intriguing. The elastic-plastic behaviour until the stage of necking can be captured well with the FE method, but beyond the onset of necking, the correspondence to reality for structure fracture may be very poor [2]. Simulations of fracture processes of structures need both accurate material models and failure criteria to address the failure elements. For steel materials under predominantly tensile loading, many researchers have proposed different failure criteria for ductile fracture, including the critical equivalent plastic strain (EPS) criterion, the Rice-Tracey-Cockcroft-Latham (RTCL) criterion [3, 4], the GTN damage criterion [5] [6] [7] [8], the Bressan-Williams-Hill (BWH) instability criterion [9], and the Forming Limit Diagram (FLD) [10, 11]. Although the BWH criterion and the RTCL damage criterion can predict the resistance versus penetration depth of steel structures under punching action satisfactorily with calibrated material models and appropriate mesh size [12], they cannot be used to reveal the influence mechanisms of material shearing behaviour on the structural responses. 

The porous plasticity-inspired micro-mechanical models are widely used in ductile fracture analysis. Some studies [5, 8, [13] [14] [15] [16] [17] [18] for ductile fracture analysis have been conducted using the popular GTN model. Gurson [5] first developed a pressure-dependent constitutive relation for yielding a porous material. Tvergaard and Needleman [8] then extended Gurson's model by improving the yield condition to incorporate the coalescence of micro-voids. Faleskog et al. [13] conducted the micromechanics calibration of void growth considering material strain hardening based on the GTN model. Kim et al. [14] investigated the influences of initial porosity and stress triaxiality on void growth and coalescence. Qian et al. [15, 16] investigated the ductile fracture behaviour of tubular joints and ductile crack extension of the API 5L X65 steel pipeline through the GTN plasticity damage model. Gholipour et al. [17] used the GTN model to simulate the ductile fracture of in situ tensile testing specimens under different stress states. Yan et al. [18] carried out numerical simulations for ductile fracture of the cold-formed S700 steel using the GTN model. 

Despite its success in simulating ductile fractures in metallic materials under medium to high stress triaxialities, the critical weakness of the GTN model lies in its inability in describing the fracture driven predominantly by shear actions under low stress triaxiality. In pure shear conditions where the stress triaxiality is zero, the void growth behaviour in the GTN model prohibits the void growth, which requires positive mean stress. To deal with this problem, some researchers have modified the GTN model according to the failure criterion of voids in the internal shear band proposed by McClintock et al. [19]. Xue [20, 21] developed an extended stress-state function to distinguish tension and shear stress states and modified the GTN model considering the shear actions of voids. Then, Nahshon and Hutchinson [22] proposed a damage law in the GTN model to account for the damage accumulation under shearing. Butcher et al. [23] numerically simulated steel straight tube hydroforming (DP600) using a modified GTN constitutive model to account for the influence of void shape and shear on coalescence. Zhou et al. [24] coupled the volumetric damage and the shear damage into the yield function and flow potential of the GTN model, which can simulate shear-induced damage by combining the damage mechanics method with the void growth model. Malcher et al. [25] assessed two different shear mechanisms proposed by Xue [20] and Nahshon and Hutchinson [22] for the prediction of crack formation. Li et al. [26] developed a modified GTN-Thomason criterion for modelling ductile failures of a metal material. He et al. [27] proposed an ISMGTN model considering shear failure and containing two independent damage mechanisms. This model was used for ductile fracture prediction of aluminium alloys under different stress states. 

The literature review shows that the modified GTN models considering the shear influences on coalescence are widely used in ductile fracture analysis under different stress states for some small-scale material specimens. Very little work for ductile failure analysis of a large-scale thin-walled structure using the above-mentioned GTN models is reported. The structures under quasi-static punch load exhibited strong shear failure behaviours in both the accidental scenarios and model tests, as mentioned in Ref. [28]. It is difficult for the current failure models to clearly investigate the influences of material shear failure characteristics on structural responses under quasi-static punch loading, such as the EPS criterion, the RTCL criterion, and the BWH instability criterion. Developing a reliable and accurate numerical failure model to simulate the shear-dominated ductile fracture of largescale, thin-walled steel structures is still a research difficulty due primarily to the challenges in determining the fracture criterion of steel materials under complex stress states. 

Therefore, this paper aims to apply the modified GTN model considering shear effects to numerical failure analysis of steel structures under quasi-static punch loading and assess the influences of the shear damage on coalescence on the resistance to penetration response of the structure. The accuracy of simulations for the ductile fracture of steel structures under quasi-static punch loading using the modified GTN model are verified by experimental results in the paper. This paper starts with the description of the modified GTN model considering shear effects proposed by Xue [20] and Nahshon and Hutchinson [22]. The modified GTN damage models presented in this section are implemented in ABAQUS software via a VUMAT subroutine. The next section presents a series of numerical simulations on the ductile fracture of in situ tensile tests under different stress states using the modified GTN model. In this section, the void-related and damage parameters in the model are calibrated based on the experimental results. The subsequent section carried out numerical simulations for fracture behaviours of non-stiffened and stiffened steel plates under quasi-static punch loading using the modified GTN model. The numerically obtained force-displacement results are compared with the experimental ones. The last section summarizes the conclusions drawn from this study. The research of this paper can provide a useful tool to simulate the failure behaviour of structures in ship collision and grounding.
section_name: Introduction

section_name: Description of the Constitutive Model

section_content: Gurson [5] initially developed a pressure-dependent yield function from an isolated spherical void in a continuum medium. Then, Tvergaard and Needleman [6] [7] [8] enhanced Gurson's model to construct the GTN model theory containing porous plasticity. In the GTN model [5] [6] [7] [8], the yield function of the GTN model is written as

where φ is the yield function, σ eq is the von Mises stress of material,

and σ m is the yield stress as a function of equivalent plastic strain in the matrix ε p m ; σ h is the macroscopic hydrostatic pressure,

q 2, q 3 are the coefficients, and q 3 = q 2 1. The equivalent void volume fraction f* is a function of the current void volume fraction f :

where f c and f f are the critical void fraction and the void volume fraction at fracture, respectively. The void coalescence will occur when the void volume fraction reaches a critical value f c, and the material will fail when the void volume fraction reaches the value of f f. Equation (2) works very well for low-stress triaxiality cases but slightly worse for high triaxiality cases [29]. 

The increase in the void volume fraction under loading comes from two parts: the original void growth and the extended void nucleation. The total change rate of the void volume fraction is written as

where. 

f growth is the change rate of voids owing to growth, and. 

f nucleation is the change rate of voids owing to the extended void nucleation. 

The matrix material is assumed to be plastically incompressible. Therefore, the void growth rate can be expressed as. 

where I is the second-order unit tensor, and. ε p is the macroscopic plastic strain rate tensor. 

A normal distribution of void nucleation versus the plastic strain is expressed as. 

where. ε p m is the plastic strain rate, and the void nucleation coefficient A controlled by strain is

where f N represents the limit value of volume fraction of all particles with potential for micro-void nucleation, ε N is the average equivalent plastic strain, s N is the standard deviation of the distribution. 

In porous materials, the plastic flow is related to the void volume fraction and the accumulated plastic strain. The matrix plastic strain is related to the macroscopic strain by the equivalence of the rate of the plastic work, which can be shown as σ :

. 

where σ is Cauchy's stress tensor. Therefore, the cumulative plastic strain can be written as. 

If the hydrostatic pressure effect is ignored, the yield function of the GTN model will be

When f = 0, the yield surface is the original Mises yield surface. The presence of a nonzero f value effectively shrinks the original Mises yield surface, as shown Figure 1. According to Equation (1), the yield surface in the σ eq -σ h plane under different equivalent void volume fractions can be obtained as shown in Figure 1. The damage accumulation forces a strong dependence of the yield surface and the plastic potential on the hydrostatic stress, in contrast to the original Mises yield surface which is independent of the mean stress. 

According to Equation (1), the yield surface in the eq σ -h σ plane under different equiv- alent void volume fractions can be obtained as shown in Figure 1. The damage accumulation forces a strong dependence of the yield surface and the plastic potential on the hydrostatic stress, in contrast to the original Mises yield surface which is independent of the mean stress.
section_name: GTN Model

section_content: The original GTN model cannot capture localization and fracture for low triaxialities in shear-dominated deformations. Under shear actions, the distortion of voids and intervoid linking promotes an effective increase in the material's internal degradation and contributes to material softening.
section_name: Modified GTN model for Shear Failure

section_content: Xue [20] further extended the original GTN model to incorporate the void shearing damage and assumed that the relative position of the void does not change, as shown in Figure 2. With the shear strain increasing, the distance between the free surface and the boundary of the representative volume element decreases. The cell in the initial undeformed configuration and in the deformed configuration is shown in Figure 2. 

In the undeformed configuration, the minimal distance between the free surface and the boundary of the cell is

The shear strain γ in the deformed configuration is

where α is the deformation angle seen in Figure 2b. The minimal distance at a deformed configuration is
section_name: Xue's Modified GTN Model

section_content: The original GTN model cannot capture localization and fracture for low triaxialities in shear-dominated deformations. Under shear actions, the distortion of voids and inter-void linking promotes an effective increase in the material's internal degradation and contributes to material softening. With the logarithm definition of strain, the artificial strain art ε is defined a In the undeformed configuration, the minimal distance between the free surface and the boundary of the cell is

The shear strain γ in the deformed configuration is

where α is the deformation angle seen in Figure 2b. The minimal distance at a deformed configuration is

With the logarithm definition of strain, the artificial strain ε art is defined as

Additionally, the failure shear strain in the shear band ε shearband for a small void volume is defined as

The damage with respect to the shearing of the void is expressed as

Using Taylor's expansion for γ = √ 3ε eq, the artificial strain is reduced to the leading terms, which can be expressed as

Therefore, Equation ( 17 ) can be rewritten as

where f = πR 2 /L 2 is the void volume fraction of the cell. Equation ( 17 ) can be only used for a 2D problem. For the 3D case, a spherical void of radius R is at the centre of the representative cell of length L. A similar relation for damage associated with void shearing can be obtained as

where f = 4πR 3 / 3L 3 is the void volume fraction for the 3D cell. The evolution of shear damage can be written using the rate form:

. D shear = q 4 f q 5 ε eq. ε eq (19) where q 4 = 1. 69 and q 5 = 1/2 for a 2D cell structure, and q 4 = 3. 72 and q 5 = 1/3 for a 3D cell structure. Xue [20]

The Lode angle related to the third invariant of the deviatoric stress tensor is an important parameter in the characterization of the stress rate effect on the ductile failure of the material. The function of the Lode angle proposed by Xue [20] is defined as

where g 0 represents the Lode angle function, and θ is the normalized Lode angle which assumes values between -1 ≤ θ ≤ 1. The end conditions of g 0 = 0 for the generalized tension (θ L = -30 • ) and g 0 = 1 for the generalized shear (θ L = 0 • ) should be satisfied, where the Lode angle θ L is

In Equation (22), χ = (s 2s 3 )/(s 1s 3 ), and s 1, s 2, and s 3 are the principal values of the deviatoric stress tensor.
section_name: Modified GTN model for Shear Failure

section_content: Nahshon and Hutchinson [22] developed a function to quantify the shearing damage effects of voids at low triaxiality, which was directly added to the total change rate of the void volume fraction, although it cannot represent a physical accumulation of the porosity. The extended shear failure term assuming linear relations of the porosity and the deviatoric plastic strain rate is written as. D shear = k s f ω s :

. 

The void volume fraction rate is expressed by the sum of the growth of voids and shear effects by. 

where k s is a defined damage parameter used for scaling the magnitude of the damage growth rate in shear loading, and ω is the function of the stress states, which is defined by

where ξ is the Lode parameter defined as

where the q is the equivalent stress, q = σ eq, and the third invariant of stress J 3 is

In Equation ( 14 ), s 1, s 2, s 3 are the principal values of the deviatoric stress tensor. A Lode angle dependence function was proposed by Nahshon and Hutchinson [22] :

where ξ is the normalized third invariant, -1 ≤ ξ ≤ 1. Figure 3 compares the evolution of the Lode angle functions proposed by Nahshon and Hutchinson [22] and Xue [20]. 

A Lode angle dependence function was proposed by Nahshon and Hutchin The shear damage evolutions can be rephrased for arbitrary loading conditi g q f D g k

The modified GTN model can simulate the fracture mode dominated by t mechanism under low stress triaxiality and account for the effective damage accum caused by shear mechanisms such as rotation and torsion of micro-voids under plex stress stage. This model can provide some physical and mechanical explana the ductile fracture of metal to a certain extent. However, there are some issues t to be studied and discussed. The parameters and damage variables in this mode be identified and calibrated through extensive numerical simulations, and the m sitivity of the model needs to be discussed. This model is mainly used in ductile analysis under different stress states for some small-scale material specimens but used for failure analysis of a large-scale thin-walled structure. The study for the tion of this model in failure analysis of large-scale thin-walled structures needs t ducted.
section_name: Nahshon and Hutchinson's (N-H) Modified GTN Model

section_content: The finite element method is widely used in failure analysis of different structures [30] [31] [32] [33]. The present section presents the numerical simulations for qu tensile tests to study the fracture mechanism of specimens with various geometri different stress states. The modified GTN damage model is implemented in ABAQ ware through using the VUMAT subroutine. The void-related parameters and s ure parameters of the modified GTN damage model elaborated in Section 2 are c The shear damage evolutions can be rephrased for arbitrary loading conditions as. D shear = g 0 q 4 f q 5 ε eq. ε eq for Xue s model

where

The modified GTN model can simulate the fracture mode dominated by the shear mechanism under low stress triaxiality and account for the effective damage accumulation caused by shear mechanisms such as rotation and torsion of micro-voids under the complex stress stage. This model can provide some physical and mechanical explanations for the ductile fracture of metal to a certain extent. However, there are some issues that need to be studied and discussed. The parameters and damage variables in this model need to be identified and calibrated through extensive numerical simulations, and the mesh sensitivity of the model needs to be discussed. This model is mainly used in ductile fracture analysis under different stress states for some small-scale material specimens but is rarely used for failure analysis of a large-scale thin-walled structure. The study for the application of this model in failure analysis of large-scale thin-walled structures needs to be conducted.
section_name: Quasi-Static Tension Tests Correlations

section_content: The finite element method is widely used in failure analysis of different types of structures [30] [31] [32] [33]. The present section presents the numerical simulations for quasi-static tensile tests to study the fracture mechanism of specimens with various geometries under different stress states. The modified GTN damage model is implemented in ABAQUS software through using the VUMAT subroutine. The void-related parameters and shear failure parameters of the modified GTN damage model elaborated in Section 2 are calibrated based on the results of tensile tests. The correlated parameters can be used for modelling the ductile fracture of a complex structure under quasi-static punch loadings in the next section.
section_name: Quasi-Static Tension Tests Correlations

section_content: In order to calibrate the parameters of the modified GTN model, seven in situ quasistatic tensile tests [32, 33] were chosen to validate the numerical results. Figure 4 shows the dimensions of seven tensile specimens with different shapes. The tensile tests for the first four specimens as shown in Figure 4a -d were conducted by Kõrgesaar et al. [32]. The tensile tests for the last two specimens as shown in Figure 4e,f were carried out by Zhang et al. [33]. The force-displacement responses as well as the failure initiation and crack propagation corresponding to each specimen were obtained in the tensile tests [32, 33]. section.
section_name: Finite Element Models for Tensile Specimens

section_content: In order to calibrate the parameters of the modified GTN model, seven in situ quasistatic tensile tests [32, 33] were chosen to validate the numerical results. Figure 4 shows the dimensions of seven tensile specimens with different shapes. The tensile tests for the first four specimens as shown in Figure 4a -d were conducted by Kõrgesaar et al. [32]. The tensile tests for the last two specimens as shown in Figure 4e,f were carried out by Zhang et al. [33]. The force-displacement responses as well as the failure initiation and crack propagation corresponding to each specimen were obtained in the tensile tests [32, 33]. [33] ; (e,f) Zhang et al. [33]. )

The engineering and true stress-strain relations of steel materials can be obtained according to the load-displacement curves from the tensile tests ST-1 and ST-2. Table 1 lists two types of material parameters of steel specimens. For material 1 [32], the true stress-strain relationship follows a Swift model:

In this equation, eq σ is the equivalent true stress, y σ is the yield stress, and k and n are the parameters of strain hardening. The strain hardening is delayed until the plastic strain reaches the plateau strain L ε. The parameter 0 ε is ( )

For material 2 [33], the equivalent true stress-strain relationship is expressed by a pow-law relation: [33] ; (e,f) Zhang et al. [33] ). 

The engineering and true stress-strain relations of steel materials can be obtained according to the load-displacement curves from the tensile tests ST-1 and ST-2. Table 1 lists two types of material parameters of steel specimens. For material 1 [32], the true stress-strain relationship follows a Swift model: In this equation, σ eq is the equivalent true stress, σ y is the yield stress, and k and n are the parameters of strain hardening. The strain hardening is delayed until the plastic strain reaches the plateau strain ε L. The parameter ε 0 is

For material 2 [33], the equivalent true stress-strain relationship is expressed by a pow-law relation:

A series of 3D models were built with the eight-node solid element with reduced integration (C3D8R) in the FE analysis, as shown in Figure 5. The refined and coarse meshes are used within and outside the central region of specimens, respectively, in order to better capture fracture behaviours of specimens and to reduce the simulation times. A series of 3D models were built with the eight-node solid element with integration (C3D8R) in the FE analysis, as shown in Figure 5. The refined and meshes are used within and outside the central region of specimens, respectively, to better capture fracture behaviours of specimens and to reduce the simulation t
section_name: Finite Element Models for Tensile Specimens

section_content: In the original GTN model, nine material parameters need to be determined ibrated based on experimental and numerical results, including coefficients of t function (q1, q2 and q q = ), void nucleation parameters (f0, fN, N ε, SN), and void fractions of coalescence and fracture (fc, ff). In order to ascertain the range of origin model material parameters, Table 2 summarizes the GTN model parameters r from a literature review. Most of the references in Table 2 adopt the same quan yield function parameters proposed by Tvergaard et al. [8], q1 = 1, q2 = 1. 5, q 0. 3 N ε =, and SN = 0. 1. For the N-H modified GTN model, the parameter ks f
section_name: Parameter Calibration of the Modified GTN Model

section_content: In the original GTN model, nine material parameters need to be determined and calibrated based on experimental and numerical results, including coefficients of the yield function (q 1, q 2 and q 3 = q 2 1 ), void nucleation parameters (f 0, f N, ε N, S N ), and void volume fractions of coalescence and fracture (f c, f f ). In order to ascertain the range of original GTN model material parameters, Table 2 summarizes the GTN model parameters retrieved from a literature review. Most of the references in Table 2 adopt the same quantities of yield function parameters proposed by Tvergaard et al. [8], q 1 = 1, q 2 = 1. 5, q 3 = 2. 25, ε N = 0. 3, and S N = 0. 1. For the N-H modified GTN model, the parameter k s for shear damage needs to be determined. For Xue's modified GTN model, the parameters q 4 and q 5 need to be calibrated. 

Table 3 lists the range of values for material parameters in the modified GTN models. Based on the numerical simulations for the tensile test (ST-2), the sensitivity analysis of the material parameters (q 1, q 2, f 0, f N, f c, f f ) of the GTN model is shown in Figure 6.
section_name: Parameter Calibration of the Modified GTN Model

section_content: Material The comparisons of numerically obtained load-displacement curves under different GTN parameters are given in Figure 6. They show that the final tension displacement increases with the increase in the parameter values (f c, f f ) or with the decrease in the parameter values (q 1, q 2, f 0, f N ). The parameter of initial void volume fraction f 0 affects the load-displacement response, as it causes a shrink in the yield surface of a material as seen in Figure 1. The other parameters (q 1, q 2, f N, f c, f f ) have no effect on the load-displacement curve prior to necking. Figure 6 also shows that the parameter of the volume fraction for void nucleation f N imposes a strong effect on the load-displacement response of specimen under uniaxial tension. Figure 7 shows the sensitivity analysis of the parameters (q 4 and q 5 ) for shear damage in Xue's modified GTN model. As the parameter q 4 decreases or the parameter q 5 increases, the final tension displacement will increase. The larger value of q 4 or the smaller value of q 5 increases the shear damage rate. D shear, which facilitates the fracture failure in the specimen. Table 3 lists the range of values for material parameters in the modified GTN models Based on the numerical simulations for the tensile test (ST-2), the sensitivity analysis o the material parameters (q1, q2, f0, fN, fc, ff) of the GTN model is shown in Figure 6. The comparisons of numerically obtained load-displacement curves under differen GTN parameters are given in Figure 6. They show that the final tension displacement in creases with the increase in the parameter values (fc, ff) or with the decrease in the param eter values (q1, q2, f0, fN). The parameter of initial void volume fraction f0 affects the loaddisplacement response, as it causes a shrink in the yield surface of a material as seen in Figure 1. The other parameters (q1, q2, fN, fc, ff) have no effect on the load-displacemen curve prior to necking. Figure 6 also shows that the parameter of the volume fraction for void nucleation fN imposes a strong effect on the load-displacement response of specimen under uniaxial tension. Figure 7 shows the sensitivity analysis of the parameters (q4 and q5) for shear damage in Xue's modified GTN model. As the parameter q4 decreases or the parameter q5 increases, the final tension displacement will increase. The larger value of q or the smaller value of q5 increases the shear damage rate Based on the sensitivity analysis of the modified GTN model parameters, Table 4 lists the calibrated material parameters which lead to the closest agreement between the numerically computed and the experimentally obtained load-displacement curves. To account for the mesh size sensitivity of the modified GTN model, FE simulations are conducted with four types of uniform element sizes in the failure zone. For the specimen ST- Based on the sensitivity analysis of the modified GTN model parameters, Table 4 lists the calibrated material parameters which lead to the closest agreement between the numerically computed and the experimentally obtained load-displacement curves. To account for the mesh size sensitivity of the modified GTN model, FE simulations are con- Based on the sensitivity analysis of the modified GTN model parameters, Table 4 lists the calibrated material parameters which lead to the closest agreement between the numerically computed and the experimentally obtained load-displacement curves. To account for the mesh size sensitivity of the modified GTN model, FE simulations are conducted with four types of uniform element sizes in the failure zone. For the specimen ST-1, the mesh lengths are chosen as 0. 3 mm, 0. 5 mm, 1 mm, and 2 mm, which correspond to four different ratios of mesh size over the plate thickness, 0. 1, 0. 167, 0. 333, and 0. 667, respectively. For the specimen ST-2, the mesh lengths are chosen as 0. 3 mm, 0. 5 mm, 1 mm, and 2 mm, which have the corresponding ratios of mesh size over the plate thickness of 0. 0952, 0. 159, 0. 317, and 0. 635.
section_name: Source

section_content: Material-1 Figure 8 shows the comparisons of engineering stress-strain curves between standard tensile tests (ST-1, ST-2) and numerical simulations with different mesh sizes based on Xue's modified GTN model. Before the peak force is reached, the simulated engineering stress-strain curves show a reasonable agreement with the experiments for all mesh sizes. A strong mesh sensitivity in the modified GTN model arises at large deformation levels. Differences between the load-displacement curves become more pronounced when the mesh size increases. The displacement at the onset of fracture with fine meshes is obviously smaller than that with medium or coarse meshes. Figure 9 compares the numerically obtained load-displacement curves based on three material models, the original GTN model, N-H modified GTN model, and Xue's modified GTN model. In Figure 9, the FE models used have a length of 1 mm. The comparisons show that the predicted final tension displacement using the original GTN model is larger than those using the other two modified GTN models, and the shear damage parameter k s of the N-H modified GTN model has little influence on the load-displacement response in the standard tension state. Figure 8 shows the comparisons of engineering stress-strain curves between standard tensile tests (ST-1, ST-2) and numerical simulations with different mesh sizes based on Xue's modified GTN model. Before the peak force is reached, the simulated engineering stress-strain curves show a reasonable agreement with the experiments for all mesh sizes. A strong mesh sensitivity in the modified GTN model arises at large deformation levels. Differences between the load-displacement curves become more pronounced when the mesh size increases. The displacement at the onset of fracture with fine meshes is obviously smaller than that with medium or coarse meshes. Figure 9 compares the numerically obtained load-displacement curves based on three material models, the original GTN model, N-H modified GTN model, and Xue's modified GTN model. In Figure 9, the FE models used have a length of 1 mm. The comparisons show that the predicted final tension displacement using the original GTN model is larger than those using the other two modified GTN models, and the shear damage parameter ks of the N-H modified GTN model has little influence on the load-displacement response in the standard tension state. For the tensile shear specimen (STS), the numerical crack initiation and propagation based on the N-H modified GTN model show a good agreement with the experimental observation. In Figure 10d, for Xue's modified GTN model, the fracture of the STS specimen is captured with a large nominal tension displacement under shear-dominated loading conditions. Compared to the N-H modified GTN model, Xue's modified GTN model behaves more conservatively, as its peak load is much larger than that of the N-H modified GTN model, and it behaves a much slower post-peak load reduction. For the N-H modified GTN model, the STS specimen fractures with a small nominal tension displacement similar to the experimental data under the same loading condition. Moreover, with the increase in the parameter ks of shear damage, the material fails quickly with a small tension displacement, as the damage rate For the tensile shear specimen (STS), the numerical crack initiation and propagation based on the N-H modified GTN model show a good agreement with the experimental observation. In Figure 10d, for Xue's modified GTN model, the fracture of the STS specimen is captured with a large nominal tension displacement under shear-dominated loading conditions. Compared to the N-H modified GTN model, Xue's modified GTN model behaves more conservatively, as its peak load is much larger than that of the N-H modified GTN model, and it behaves a much slower post-peak load reduction. For the N-H modified GTN model, the STS specimen fractures with a small nominal tension displacement similar to the experimental data under the same loading condition. Moreover, with the increase in the parameter ks of shear damage, the material fails quickly with a small tension displacement, as the damage rate For the tensile shear specimen (STS), the numerical crack initiation and propagation based on the N-H modified GTN model show a good agreement with the experimental observation. In Figure 10d, for Xue's modified GTN model, the fracture of the STS specimen is captured with a large nominal tension displacement under shear-dominated loading conditions. Compared to the N-H modified GTN model, Xue's modified GTN model behaves more conservatively, as its peak load is much larger than that of the N-H modified GTN model, and it behaves a much slower post-peak load reduction. For the N-H modified GTN model, the STS specimen fractures with a small nominal tension displacement similar to the experimental data under the same loading condition. Moreover, with the increase in the parameter k s of shear damage, the material fails quickly with a small tension displacement, as the damage rate. D shear associated with void shearing is positively correlated with the parameter k s. Figure 13 illustrates the evolution of the stress triaxiality with respect to the equivalent plastic strain for each tensile specimen. The stress triaxiality is defined as the ratio of the macroscopic hydrostatic pressure and equivalent effective stress of material,

. The evolution of the stress triaxiality is evaluated at the element located in the region (marked by the red dot in Figure 13 ) where the crack initiates. The different stress triaxiality levels for the same material are observed from the different tensile tests with different geometries. As shown in Figure 13, for specimens ST-1, ST-2, CH, and NT-1, the stress triaxialities at low plastic strain levels obtained from three GTN damage models are very close. Moreover, for the standard tensile specimens (ST-1, ST-2), notched specimens (NT-1, NT-2), and central-hole specimen (CH), the stress triaxialities exceed the value of 1/3 corresponding to the uniaxial tension. In contrast, the initial stress triaxialities in the sheared tensile specimen (STS) remain below the uniaxial tension value of 1/3 at a small deformation level. The predicted stress triaxialities versus equivalent plastic strain for all specimens based on three types of GTN models have some differences as the damage evolves in the element with increasing deformation. The damage evolution associated with the shearing of voids has significant influences on the stress triaxiality versus equivalent plastic strain, especially for the shear-tension-dominated loading conditions. Figure 13 illustrates the evolution of the stress triaxiality with respect to the equivalent plastic strain for each tensile specimen. The stress triaxiality is defined as the ratio of the macroscopic hydrostatic pressure and equivalent effective stress of material, σ h /σ eq. The evolution of the stress triaxiality is evaluated at the element located in the region (marked by the red dot in Figure 13 ) where the crack initiates. The different stress triaxiality levels for the same material are observed from the different tensile tests with different geometries. As shown in Figure 13, for specimens ST-1, ST-2, CH, and NT-1, the stress triaxialities at low plastic strain levels obtained from three GTN damage models are very close. Moreover, for the standard tensile specimens (ST-1, ST-2), notched specimens (NT-1, NT-2), and centralhole specimen (CH), the stress triaxialities exceed the value of 1/3 corresponding to the uniaxial tension. In contrast, the initial stress triaxialities in the sheared tensile specimen (STS) remain below the uniaxial tension value of 1/3 at a small deformation level. The predicted stress triaxialities versus equivalent plastic strain for all specimens based on three types of GTN models have some differences as the damage evolves in the element with increasing deformation. The damage evolution associated with the shearing of voids has significant influences on the stress triaxiality versus equivalent plastic strain, especially for the shear-tension-dominated loading conditions. Figure 13 illustrates the evolution of the stress triaxiality with respect to the equivalent plastic strain for each tensile specimen. The stress triaxiality is defined as the ratio of the macroscopic hydrostatic pressure and equivalent effective stress of material,

. The evolution of the stress triaxiality is evaluated at the element located in the region (marked by the red dot in Figure 13 ) where the crack initiates. The different stress triaxiality levels for the same material are observed from the different tensile tests with different geometries. As shown in Figure 13, for specimens ST-1, ST-2, CH, and NT-1, the stress triaxialities at low plastic strain levels obtained from three GTN damage models are very close. Moreover, for the standard tensile specimens (ST-1, ST-2), notched specimens (NT-1, NT-2), and central-hole specimen (CH), the stress triaxialities exceed the value of 1/3 corresponding to the uniaxial tension. In contrast, the initial stress triaxialities in the sheared tensile specimen (STS) remain below the uniaxial tension value of 1/3 at a small deformation level. The predicted stress triaxialities versus equivalent plastic strain for all specimens based on three types of GTN models have some differences as the damage evolves in the element with increasing deformation. The damage evolution associated with the shearing of voids has significant influences on the stress triaxiality versus equivalent plastic strain, especially for the shear-tension-dominated loading conditions.
section_name: Material Type

section_name: Steel-Plated Structures under Quasi-Static Punch Loading

section_content: Zhang et al. [33] presented the experimental study on non-stiffened plates under quasi-static punch loading. Figure 14 illustrates the experimental setup and dimensions of indenter and plate. The thickness of the plate is equal to 3. 15 mm. In the test, each of the four edges of the plate specimen had clamp supports between two thick steel plates, which were fixed with the strong base. The plate specimen was loaded at a rate of about 10 mm/min at the centre by hydraulic cylinders. The material used for the non-stiffened plate remained the same as Material-2 in Table 1.
section_name: Simulations of Non-Stiffened Plate Penetration Test

section_content: Zhang et al. [33] presented the experimental study on non-stiffened plates under quasi-static punch loading. Figure 14 illustrates the experimental setup and dimensions of indenter and plate. The thickness of the plate is equal to 3. 15 mm. In the test, each of the four edges of the plate specimen had clamp supports between two thick steel plates, which were fixed with the strong base. The plate specimen was loaded at a rate of about 10 mm/min at the centre by hydraulic cylinders. The material used for the non-stiffened plate remained the same as Material-2 in Table 1. 

The non-stiffened plate penetration test was analysed using the ABAQUS software with the calibrated modified GTN model given in previous sections. The material parameters of the plate specimen (material-2) are presented in Table 1, and the calibrated input parameters for GTN models are given in Table 4. The eight-node solid element with reduced integration (C3D8R) was used in modelling the plate specimen, as shown in Figure 15. The plate specimen in the contact region was meshed with three different characteristic solid element sizes, 1 mm, 3 mm, and 5 mm. All degrees of freedoms for nodes at the boundary edges of the plate were constrained in the numerical model to simulate the fully clamped boundary condition in the penetration test, without modelling the detailed clamping setup. The indenter was modelled with rigid material because no deformations were observed in the test. The surface-to-surface contact (explicit) defines the interaction between the indenter and plate specimen. The friction coefficient in the contact interaction is defined as 0. 2 [33]. The part of the plate specimen in the contact region where the crack initiates and propagates is simulated with fine solid elements. The non-stiffened plate penetration test was analysed using the ABAQUS sof with the calibrated modified GTN model given in previous sections. The material p eters of the plate specimen (material-2) are presented in Table 1, and the calibrated parameters for GTN models are given in Table 4. The eight-node solid element wi duced integration (C3D8R) was used in modelling the plate specimen, as shown in F 15. The plate specimen in the contact region was meshed with three different charact solid element sizes, 1 mm, 3 mm, and 5 mm. All degrees of freedoms for nodes boundary edges of the plate were constrained in the numerical model to simulate the clamped boundary condition in the penetration test, without modelling the de clamping setup. The indenter was modelled with rigid material because no deform were observed in the test. The surface-to-surface contact (explicit) defines the inter between the indenter and plate specimen. The friction coefficient in the contact inter is defined as 0. 2 [33]. The part of the plate specimen in the contact region where the initiates and propagates is simulated with fine solid elements. The non-stiffened plate penetration test was analysed using the ABAQUS sof with the calibrated modified GTN model given in previous sections. The material p eters of the plate specimen (material-2) are presented in Table 1, and the calibrated parameters for GTN models are given in Table 4. The eight-node solid element wi duced integration (C3D8R) was used in modelling the plate specimen, as shown in F 15. The plate specimen in the contact region was meshed with three different charact solid element sizes, 1 mm, 3 mm, and 5 mm. All degrees of freedoms for nodes boundary edges of the plate were constrained in the numerical model to simulate the clamped boundary condition in the penetration test, without modelling the de clamping setup. The indenter was modelled with rigid material because no deform were observed in the test. The surface-to-surface contact (explicit) defines the inter between the indenter and plate specimen. The friction coefficient in the contact inter is defined as 0. 2 [33]. The part of the plate specimen in the contact region where the initiates and propagates is simulated with fine solid elements. Figure 16 shows experimental load-displacement curve and numerical load-displacement curves computed from different meshes, GTN models, and damage parameters of k s. In the test, the contact between the plate specimen and the indenter has three stages as the indenter penetrates through the steel plate: (a) hemisphere-plate contact, (b) cylinder-plate contact, and (c) cone-plate contact, which correspond to different shapes of the three segments of the indenter shown in Figure 14. The whole penetration processes of the nonstiffened plate specimen are categorized into two fracture processes: the crack initiation process and crack propagation process. In the first contact stage (hemisphere-plate contact stage), the experimentally measured load decreases sharply when the crack initiates and propagates with the increase in the penetration depth. In the next contact stage (cylinderplate contact stage), the measured load remains at a very low level with the increase in penetration depth. In the last contact stage, when the cone part of the indenter punches the ruptured plate, the load increases remarkably, as shown in Figure 16. Figure 16a and Table 5 compare the numerical load-displacement curves obtained from the original GTN model, N-H modified GTN model, and Xue's modified GTN model. The comparisons show that the numerically computed load-displacement response based on the N-H modified GTN model agrees most closely with the experimental data for all three stages of contacts. The modified GTN models considering shear effect show better predictions on the peak load and post-peak response in the experiment than the original GTN model. Figure 16b, c present the numerical load-displacement curves with different mesh sizes obtained from the N-H modified GTN model and Xue's modified GTN model, respectively. The plate specimen in the contact region is meshed with three different characteristic element sizes, 1 mm, 3 mm, and 5 mm. The mesh sizes of plate model have a marginal effect on the numerical load-displacement relations. For Xue's modified GTN model, the peak force simulated with coarse mesh size (5 mm) is closest to the experimental value among the three meshes, whereas for the N-H modified GTN model, the fine mesh size (1 mm and 3 mm) has best estimation on the experimental response. Figure 16d elucidates the influences of the damage parameter k s in the N-H modified GTN model, which shows a negligible effect on the load-displacement responses. Figure 17 compares the experimentally and numerically obtained fracture shapes of the plate specimen. The N-H modified GTN damage model with k s = 4 was used in the simulation. In the penetration test, the plate specimen experienced plastic deformations, necking, crack initiation, and propagation, which are also observed in the numerical simulation of the N-H modified GTN damage model. At first, necking appears on the plate with the increase in plastic deformations in the contact region, and the plate subsequently experiences fracture accompanied by a rapid load decrease. Several cracks subsequently emerge in the plate around the necking circle with the punch of the indenter, as shown in Figure 17b. The ruptured specimen separates into several petal pieces as the crack propagates caused by the punching of the conical part of the indenter bar as shown in Figure 17d. Figure 18 shows the fracture shapes of the damaged plate specimen obtained from numerical simulations with different GTN models and mesh sizes. Due to the different failure criteria, the fracture shapes predicted by the three GTN models are quite different. The comparisons in Figures 17 and 18 show that the N-H modified GTN model with a fine mesh size has a most reasonable prediction of crack initiation and propagation in a plate under punch loading.
section_name: Steel-Plated Structures under Quasi-Static Punch Loading 4.1. Simulations of Non-Stiffened Plate Penetration Test

section_content: Kõrgesaar et al. [32] reported a quasi-static penetration test for a laser-welded stiffened plate performed in the Strength of Material laboratory at Aalto University. Figure 19 shows the experimental setup and fracture shapes of stiffened plate specimen under punch loading. The flat-bar stiffeners with a 30 mm web height and 120 mm spacing were welded on the panel by the laser-welding technique. The in-plane dimension of the stiffened panel was 1200 × 1200 mm, and the thickness of plate and stiffeners equalled 3 mm, as illustrated in Figure 20. The edge of the stiffened plate specimen is fully clamped between two supporting thick plates made from steel S355. A hemispherical indenter with a
section_name: Simulations of Laser-Welded Stiffened Plate Penetration Test

section_content: Kõrgesaar et al. [32] reported a quasi-static penetration test for a laser-welded stiffened plate performed in the Strength of Material laboratory at Aalto University. Figure 19 shows the experimental setup and fracture shapes of stiffened plate specimen under punch loading. The flat-bar stiffeners with a 30 mm web height and 120 mm spacing were welded on the panel by the laser-welding technique. The in-plane dimension of the stiffened panel was 1200 × 1200 mm, and the thickness of plate and stiffeners equalled 3 mm, as illustrated in Figure 20. The edge of the stiffened plate specimen is fully clamped between two supporting thick plates made from steel S355. A hemispherical indenter with a radius of 120 mm was loaded to punch through the centre of the stiffened plate. The material used for the stiffened plate remained the same as material-1 in Table 1. 

The simulation method for a stiffened plate penetration test is the same with that of a non-stiffened plate penetration test. The calibrated modified GTN model is implemented in ABAQUS via a user-defined VUMAT material subroutine. Tables 1 and 4 list the material parameters of the stiffened plate specimen (material-1) and the calibrated input parameters for GTN models, respectively. The eight-node solid element with reduced integration (C3D8R) was employed in modelling the stiffened plate specimen. The stiffened plate specimen in the contact region was meshed with four types of solid element sizes, 3 mm, 4. 5 mm, 6 mm, and 12 mm. All degrees of freedom for nodes on the clamped plate edges remained constrained in the numerical model. The surface-to-surface contact defines the interaction between the indenter and plate specimen. The friction coefficient in the contact interaction is defined as 0. 23 [32]. The part of panels and stiffeners in the contact region where the crack initiates and propagates were simulated with fine solid elements, as shown in Figure 21. 

radius of 120 mm was loaded to punch through the centre of the stiffened plate. The material used for the stiffened plate remained the same as material-1 in Table 1. The simulation method for a stiffened plate penetration test is the same with that of a non-stiffened plate penetration test. The calibrated modified GTN model is implemented in ABAQUS via a user-defined VUMAT material subroutine. Tables 1 and 4 list the material parameters of the stiffened plate specimen (material-1) and the calibrated input parameters for GTN models, respectively. The eight-node solid element with reduced integration (C3D8R) was employed in modelling the stiffened plate specimen. The stiffened plate specimen in the contact region was meshed with four types of solid element sizes, 3 mm, 4. 5mm, 6mm, and 12 mm. All degrees of freedom for nodes on the clamped plate edges remained constrained in the numerical model. The surface-to-surface contact defines the interaction between the indenter and plate specimen. The friction coefficient in radius of 120 mm was loaded to punch through the centre of the stiffened plate. The material used for the stiffened plate remained the same as material-1 in Table 1. The simulation method for a stiffened plate penetration test is the same with that of a non-stiffened plate penetration test. The calibrated modified GTN model is implemented in ABAQUS via a user-defined VUMAT material subroutine. Tables 1 and 4 list the material parameters of the stiffened plate specimen (material-1) and the calibrated input parameters for GTN models, respectively. The eight-node solid element with reduced integration (C3D8R) was employed in modelling the stiffened plate specimen. The stiffened plate specimen in the contact region was meshed with four types of solid element sizes, 3 mm, 4. 5mm, 6mm, and 12 mm. All degrees of freedom for nodes on the clamped plate edges remained constrained in the numerical model. The surface-to-surface contact defines the interaction between the indenter and plate specimen. The friction coefficient in Figure 22 illustrates the numerically predicted force-displacement curves with different GTN damage models and mesh sizes and compares them with experiment results. The comparisons of different GTN models in Figure 22a and Table 6 show that the numerically predicted results based on the modified GTN models considering shear failure are closer to experimental results than those of the original GTN model. Figure 22b,c illustrate the mesh size sensitivity of modified GTN damage models. The peak load and the final displacement decrease slightly with the decreasing mesh size of a stiffened plate model. Xue's and the N-H modified GTN models with a mesh size of 3 mm have agreed most closely with the experimentally measured load-deformation curves among the four element sizes considered. The mesh sensitivity analysis shows the coarse mesh elements may lead to the over-prediction of load resistance and fail to capture the strain concentrations and hence the crack initiation. Figure 22d the contact interaction is defined as 0. 23 [32]. The part of panels and stiffeners in the c tact region where the crack initiates and propagates were simulated with fine solid ments, as shown in Figure 21. Figure 22 illustrates the numerically predicted force-displacement curves with ferent GTN damage models and mesh sizes and compares them with experiment resu The comparisons of different GTN models in Figure 22a and Table 6 show that the num ically predicted results based on the modified GTN models considering shear failure closer to experimental results than those of the original GTN model. Figure 22b,c illust the mesh size sensitivity of modified GTN damage models. The peak load and the fi displacement decrease slightly with the decreasing mesh size of a stiffened plate mo Xue's and the N-H modified GTN models with a mesh size of 3 mm have agreed m closely with the experimentally measured load-deformation curves among the four ment sizes considered. The mesh sensitivity analysis shows the coarse mesh elements m lead to the over-prediction of load resistance and fail to capture the strain concentrati and hence the crack initiation. Figure 22d tact region where the crack initiates and propagates were simulated with fine solid ele-ments, as shown in Figure 21. Figure 22 illustrates the numerically predicted force-displacement curves with different GTN damage models and mesh sizes and compares them with experiment results The comparisons of different GTN models in Figure 22a and Table 6 show that the numerically predicted results based on the modified GTN models considering shear failure are closer to experimental results than those of the original GTN model. Figure 22b, c illustrate the mesh size sensitivity of modified GTN damage models. The peak load and the final displacement decrease slightly with the decreasing mesh size of a stiffened plate model Xue's and the N-H modified GTN models with a mesh size of 3 mm have agreed most closely with the experimentally measured load-deformation curves among the four element sizes considered. The mesh sensitivity analysis shows the coarse mesh elements may lead to the over-prediction of load resistance and fail to capture the strain concentrations and hence the crack initiation. Figure 22d Figure 23 presents the numerically predicted post-fracture shape of the stiffened plate with different GTN models and mesh sizes and some test photographs from Kõrgesaar et al. [33]. In the test, cracks firstly initiated on the panel and propagated parallel to the central stiffener, followed by the crack initiation and propagation on the panel in the transverse direction across the stiffener, as shown in Figure 23(a1-a3 ). Figure 23(b1-b3 ) compares the numerically computed post-fracture shape of stiffened plate from three types of GTN models. No crack occurs on the stiffened panels based on the original GTN model, whereas some small circular cracks emerge around the central stiffener in the contact region for Xue's and the N-H modified GTN models. The comparisons in Figure 23(b1-b3 ) show that the shear-driven failure is a key mechanism in the ductile fracture of plates under quasi-static punch loading. The GTN damage model requires a damage criterion for the shearing of voids to accurately reflect the failure modes in these plates. For Xue's modified GTN model, the mesh sizes of the stiffened plate model have some influence on the fracture shape, and the finer mesh facilitates the fracture in the plate structures, as shown in Figure 23(c1-c3 ). Figure 24 shows the cracking path of the stiffened plate model captured from numerical simulations based on modified GTN models. For Xue's modified GTN model, the crack first initiates and propagates on the panels in the transverse direction perpendicular to the stiffener, followed by crack initiation and propagation in the longitudinal direction along the stiffener, as shown in Figure 24a. For the N-H modified GTN model, the crack first initiates on the panels and propagates paralle to the stiffener, with subsequent crack initiation and propagation in the transverse direc- Figure 23 presents the numerically predicted post-fracture shape of the stiffened plate with different GTN models and mesh sizes and some test photographs from Kõrgesaar et al. [33]. In the test, cracks firstly initiated on the panel and propagated parallel to the central stiffener, followed by the crack initiation and propagation on the panel in the transverse direction across the stiffener, as shown in Figure 23(a1-a3 ). Figure 23(b1-b3 ) compares the numerically computed post-fracture shape of stiffened plate from three types of GTN models. No crack occurs on the stiffened panels based on the original GTN model, whereas some small circular cracks emerge around the central stiffener in the contact region for Xue's and the N-H modified GTN models. The comparisons in Figure 23(b1-b3 ) show that the shear-driven failure is a key mechanism in the ductile fracture of plates under quasi-static punch loading. The GTN damage model requires a damage criterion for the shearing of voids to accurately reflect the failure modes in these plates. For Xue's modified GTN model, the mesh sizes of the stiffened plate model have some influence on the fracture shape, and the finer mesh facilitates the fracture in the plate structures, as shown in Figure 23 (c1-c3). Figure 24 shows the cracking path of the stiffened plate model captured from numerical simulations based on modified GTN models. For Xue's modified GTN model, the crack first initiates and propagates on the panels in the transverse direction perpendicular to the stiffener, followed by crack initiation and propagation in the longitudinal direction along the stiffener, as shown in Figure 24a. For the N-H modified GTN model, the crack first initiates on the panels and propagates parallel to the stiffener, with subsequent crack initiation and propagation in the transverse direction across the stiffener. The stiffener finally ruptures completely as shown in Figure 24b
section_name: Simulations of Laser-Welded Stiffened Plate Penetration Test

section_content: This paper aims to identify an appropriate modified GTN model to simulate the shear-driven fracture failure in steel plates under quasi-static punch loading. Consequently, the original GTN model, Xue's modified GTN model, and the N-H modified model are investigated and compared to quantify the effect of shear failure on the ductile fracture and the load-displacement responses. The material and damage parameters in the different GTN models are calibrated through the laboratory material tests using specimens with varying levels of shearing and stress triaxiality. The sensitivity of a number of important parameters including the material damage parameters and the element sizes is also examined. The calibrated damage parameters for the various GTM models by the laboratory material tests are subsequently used to investigate the large-scale steel plates under punch loading. Some major concluding remarks of this paper are drawn as follows:

(1) The sensitivity analysis of void-related parameters (q 1, q 2, f 0, f N, f c, f f ) and shear damage parameters (q 4, q 5, k s ) indicate different levels of effects on the force-displacement responses for a standard tensile laboratory material test. The parameters (q 1, q 2, f N, f c, f f ) have negligible effects on the force-displacement curves before necking occurs, whereas the parameter of the initial void volume fraction f 0 affects the forcedisplacement relation before necking because a nonzero f 0 leads to a shrink in the material yield surface. The parameters of q 1, q 2, f 0, f N, and q 4 have a negative correlation and the parameters of f c, f f, and q 5 have a positive correlation with the tension displacement at the onset of fracture. The shear damage parameter k s of the N-H modified GTN model has a very minor effect on the force-displacement response in uniaxial tension. (2) The numerical prediction based on Xue's or the N-H modified GTN models provides closer estimations on the crack initiation and the load-displacement relationships observed from the tensile tests with various geometries and stress triaxialities and shearing conditions. The damage evolution associated with the shearing of voids has strong influences on the stress triaxiality versus equivalent plastic strain under complex stress states, especially for the shear-driven fracture failure. A strong mesh sensitivity to the failure criterion of the modified GTN model exists in the necking and fracture initiation stages of tests. The differences between the load-displacement curves become more pronounced when the mesh size increases. (3) The modified GTN model considering the shear effect with calibrated parameters is applied to simulate steel plate structures under quasi-static punch loading. The experimental results captured from non-stiffened plate and stiffened plate penetration tests are compared well with numerical results based on the N-H modified GTN models. The N-H modified model with k s = 4 provides much better prediction to the experimental observations compared with the original GTN model and Xue's modified model. The results demonstrated that the evolution law for the shear failure of voids in the modified GTN model is very important for simulating the ductile fracture of structures under punch loading. 

The good agreement between the experimental and numerical results shows that the modified GTN model considering the shear failure can accurately simulate cracking evolution and fracture behaviour in large-scale, thin-walled steel structures. The modified GTN model can be further applied for numerical analysis of failure behaviours of structures in ship collision and grounding events. The study of this paper can provide an effective and accurate numerical calculation procedure for the crashworthiness of a structure under localised load. But for the modified GTN models, it is still a difficult issue to give some quantitative expressions or explanations to calculate the void-related and shear damage parameters, such as f 0, f N, fc, f f, and so on, which need in-depth research in the future.
section_name: Concluding Remarks

section_content: Funding: This research was funded by the program of China Scholarship Council (CSC) (No. 201906950096 ), Foshan Xianhu Laboratory of the Advanced Energy Science and Technology Guangdong Laboratory ( XHT2020-002 ), the National Natural Science Foundation of China (no. 52273080 ), and the Fundamental Research Funds for the Central Universities (WUT: 2022IVA175 and WUT: 2023IV018h ). 

Institutional Review Board Statement: Not applicable. 

Informed Consent Statement: Not applicable.

section_name: Data Availability Statement:

section_content: The data presented in this study are available on request from the corresponding author.
section_name: Data Availability Statement:

section_name: Conflicts of Interest:

section_content: The authors declare no conflict of interest.
section_name: Conflicts of Interest: