# Comparative Study on Various Ductile Fracture Models for Marine Structural Steel EH36

## Article information

## Abstract

It is important to obtain reasonable predictions of the extent of the damage during maritime accidents such as ship collisions and groundings. Many fracture models based on different mechanical backgrounds have been proposed and can be used to estimate the extent of damage involving ductile fracture. The goal of this study was to compare the damage extents provided by some selected fracture models. Instead of performing a new series of material constant calibration tests, the fracture test results for the ship building steel EH36 obtained by Park et al. (2019) were used which included specimens with different geometries such as central hole, pure shear, and notched tensile specimens. The test results were compared with seven ductile fracture surfaces: Johnson-Cook, Cockcroft-Latham-Oh, Bai-Wierzbicki, Modified Mohr-Coulomb, Lou-Huh, Maximum shear stress, and Hosford-Coulomb. The linear damage accumulation law was applied to consider the effect of the loading path on each fracture surface. The Swift-Voce combined constitutive model was used to accurately define the flow stress in a large strain region. The reliability of these simulations was verified by the good agreement between the axial tension force elongation relations captured from the tests and simulations without fracture assignment. The material constants corresponding to each fracture surface were calibrated using an optimization technique with the minimized object function of the residual sum of errors between the simulated and predicted stress triaxiality and load angle parameter values to fracture initiation. The reliabilities of the calibrated material constants of B-W, MMC, L-H, and HC were the best, whereas there was a high residual sum of errors in the case of the MMS, C-L-O, and J-C models. The most accurate fracture predictions for the fracture specimens were made by the B-W, MMC, L-H, and HC models.

## 1. Introduction

Materials can be divided into ductile materials and brittle materials according to plastic behavior and characteristics. The brittle materials hardly exhibit plastic deformation, while the ductile materials exhibit plastic deformation accompanied by strain hardening and necking phenomenon. Mild and high strength steels mainly used in ships and marine structures and aluminum alloys used in small and medium sized fishing vessels and marine leisure vessels are involved in the category of ductile materials. The information on the fracture behavior of ductile materials is essential to minimize the accidental damages such as collision, stranding and explosion. Researchers in the shipbuilding and marine industry, as well as the automotive industry and the aerospace industry, have proposed fracture models made of various ductile materials in order to clarify the fracture behavior of the ductile materials.

The model mainly used for fracture prediction of ductile materials is divided into three types according to the mechanical viewpoint. The first is Gurson model-based fracture model. McClintock (1968) and Rice and Tracey (1969) proposed the void growth based damaged model that defined the fracture behavior of ductile materials as a series of processes of nucleation, growth and coalescence of micro void in materials and structures from a micro prospective. Gurson (1977) introduced void volume fraction as a damage parameter and proposed porous plasticity fracture model. Hancock and Brown (1983) have demonstrated the theory of micro void growth through the experiments on round bar notch specimens. Tvergaard and Needleman (1984) proposed a Gurson-Tvergaard-Needleman yield function and a fracture model (GTN model) to control the growth rate of micro void in the Gurson model. After that, many researchers have proposed fracture models based on the Gueson model, and these fracture models are used to numerically represent the damage process of materials. Second, continuum damage mechanics (CDM) has a similar concept to the GTN model, but based on the basic theory of continuum mechanics, the damage inside material is expressed as an internal parameter from a macro perspective (Lemaitre, 1985). The CDM model implements the damage inside material based on the increase of the plastic strain through the change of material stiffness (i.e., the reduction of the flow stress). Therefore, the plastic strain for damage initiation and damage evolution pattern are main material constants. The damage initiation and evolution could be represented only as function of plastic strains in the CDM. Recently, however, the CDM is being developed to represent damage as functions like characteristic displacements.

Finally, in the phenomenological model, fracture is considered to occur when the damage indicator, expressed as stress or strain at a point in the material, reaches threshold value. Unlike the two models mentioned above, in the phenomenological model, damage does not affect the material stiffness, which make it easier to calibrate material constants. Therefore, recently, many researchers are applying phenomenological models to the fracture prediction of ductile materials.

Based on pososity theory, Rice and Tracey (1969) expressed the damage as a function of the stress triaxiality, which is defined as the ratio of hydrostatic stress and von Mises equivalent stress. LeRoy et al. (1981) observed the growth of void and expressed it as a function of hydrostatic stress and maximum principal stress. Cockcroft and Latham (1968) proposed Cockcroft-Latham (C-L) expressed as the accumulation of the maximum principal stress to the increment of equivalent plastic strain, and Oh et al. (1979) modified the C-L model and expressed as a ratio of von Mises equivalent stress and maximum principal stress. Clift et al. (1990) proposed fracture model expressed as the accumulation of equivalent stress. Bao and Wierzbicki (2004) proved that the existing phenomenological model does not approximate the fracture strain prediction of experimental data through experiment at wide range of stress triaxiality condition including low stress triaxiality such as pure shear and compression as well as high stress triaxiality which has been mainly applied. After that, multiple researchers (Xue, 2007; Bai and Wierzbicki, 2008) introduced, as a variant of fracture model, Lode angle expressed as third deviatoric stress invariant along with stress triaxiality. Bai and Wierzbicki (2008) proposed second-degree polynomial with stress triaxiality and Lode angle as functions. Bai and Wierzbicki (2010) proposed the Modified Mohr-Coulomb (MMC) model expressed with stress triaxiality and Lode angle by deriving the Mohr-Coulomb model as the combination of vertical stress and shear stress acting on fracture surface. In similar way, Mohr and Marcadet (2015) proposed Hosford-Coulomb (HC) model expressed with Mohr-Coulomb fracture condition and Hosford equivalent stress. Lou et al. (2012) proposed Lou-Huh (L-H) fracture model based the growth of micro void as physical meanings, and Park et al. (2015) modified the L-H fracture model by applying Hill's 47 yield function to have anisotropic influence of a plate included in the model.

This paper will evaluate various phenomenological fracture models using some of the precedent studies (Park et al., 2019). To this end, the theoretical background of the fracture model is reviewed, and each fracture model is expressed as fracture strain with functions of stress triaxiality and Lode angle parameter.

## 2. Theoretical background

### 2.1 Stress state parameters

The stress state of isotropic materials can be expressed with stress triaxiality (*η*) and Lode angle parameter (*I*_{1}) of stress tensor (*σ*), and second invariant (*J*_{2}) and third invariant (*J*_{3}) of deviatoric stress tensor (s) (Eqs. (3)-(6)). The range of the Lode angle parameter is *η*=1/3, *η*=2/3, *η*⩽2/3).

Fig. 2 shows a cylindrical coordinate system based on the center of a specific deviatoric stress plane on the principal stress coordinate system (*σ*_{1}, *σ*_{2}, *σ*_{3}). The stress on the deviatoric stress plane perpendicular to the yield potential can be expressed as von Mises equivalent stress (*σ _{m}*) (Eq. (1)). The principal stress can be expressed with von Mises equivalent stress, stress triaxiality, and Lode angle parameter on the geometrical diagram as shown in Eqs. (8)-(10) (Mohr and Marcadet, 2015).

*f*_{1}, *f*_{2}, *f*_{3} are expressed as the functions of the Lode angle parameter as shown in Eq. (11) to Eq. (13).

### 2.2. Fracture model

In this study, total seven models (maximum shear stress, Cockcroft-Latham-Oh, Johnson-Cook, Bai-Weirzbicki, modified Mohr-Coulomb, Lou-Huh and Hosford-Coulomb) are examined, and stress strain (

#### 2.2.1 Maximum shear stress (MSS) fracture model

The maximum shear stress yield criterion are known to simulate the plastic deformation of metals accurately and have been used for a long time in various industries such as construction and civil engineering. According to the maximum shear stress yield criterion, yield is considered to occur when the maximum shear stress (*τ _{max}*) reaches a threshold value (

*τ*, yield shear strength) as shown in Eq. (14). Bai and Wierzbicki (2010) defined the stress strain with stress triaxiality and Lode angle parameters as functions as shown in Eq. (15) by deriving the maximum shear stress yield criterion based on strain. The material constants of the MSS model are composed of

_{s}*A*,

*n*,

*τ*. Here,

_{s}*A*and

*n*are the material constants of the constitutive equation (Eq. (36)) representing the flow stress. Therefore, the material constant of the MSS model is only

*τ*.

_{s}#### 2.2.2 Cockcroft-Latham-Oh (C-L-O) fracture model

Cockcroft and Latham (1968) proposed the Cockcroft-Latham (C-L) model that defines the accumulation of the maximum principal stress (*σ*_{1}) based on the increment of equivalent stress strain as fracture criterion (Eq. (16)). The C-L model is mainly used to predict the fracture occurrence generated during the forming process of sheet ductile materials. Oh et al. (1979) defined Cockcroft-Latham-Oh (C-L-O) model that the C-L model was non-dimensionized into the von Mises equivalent stress (Eq. (17)). In Eqs. (16)-(17), *C*_{1}, as strain energy and non-dimensional strain energy, is the material constant of C-L and C-L-O models. In case of proportional loading condition, non-dimensional energy *C*_{1} is not the accumulated value in Eq. (18). In addition, Eq. (8) can be substituted into Eq. (18) to derive Eq. (19), which is the final form of the C-L-O model. In Eqs. (16)-(19), Symbol < > means to take 0 if the value is negative.

#### 2.2.3 Johnson-Cook (J-C) fracture model

Johnson and Cook (1985) proposed the Johnson-Cook (J-C) fracture model composed of three material constants (*c*_{0}, *c*_{1}, *c*_{2}) and stress triaxiality (Eq. (20)). Even though the J-C fracture model does not include Lode angle effect, it has been used by many researchers due to convenient calibration of material constants. In the shipbuilding industry, Tornqvist (2003) derived the material constants of the JC fracture model through tensile tests on marine steel, and Choung et al. (2011) presented the J-C model for high stress triaxiality through the experiment on notched tension specimen of EH36 (steel grade). Min and Cho (2012) applied the J-C model to numerical analysis of drop object fracture experiment and derived material constants.

#### 2.2.4 Bai-Wierzbicki (B-W) fracture model

Bai and Wierbicki (2008) defined the fracture strain with quadratic functions of stress triaxiality and Lode angle parameter in proportional load as shown in Eq. (21). The fracture model is expressed with the material constants resulted from tensile (*D*_{1}, *D*_{2}), shear (*D*_{3}, *D*_{4}) and compression (*D*_{5}, *D*_{6}), respectively. The experiment in compression area was not included in the precedent research. Therefore, Eq. (21) can be expressed as Eq. (22), assuming that the tensile and compression terms are symmetry each other. Park et al. (2018) proposed the material constants of the B-W model through fracture experiments on steel grade, and verified the fracture prediction accuracy of the B-W model through fracture experiments and small structural tests using non-stiffening plate.

#### 2.2.5 Lou-Huh (L-H) fracture model

Lou et al. (2012) proposed the empirical model based on each process of nucleation, growth and coalescence of micro void as physical meanings from a micro-prospective on mild fracture. The fracture model is composed of terms to represent nucleation, growth and coalescence of void. In the equation for this model, the void nucleation is expressed on the basis of equivalent plastic strain. The L-H model represented the hydrostatic stress triaxiality normalized as von Mises equivalent stress, and considered Bao and Wierzbicki's research result (2005) that fracture did not occur below *η*=-1/3. The coalescence of micro-void was simulated by using the maximum shear stress and normalized to von Mises equivalent stress. The L-H model is expressed as shown in Eq. (22) and composed of three model constants (*C*_{1}, *C*_{2}, *C*_{3}). The L-H model is mainly used for the molding process of high tensile steel of thin plate and the evaluation of vehicle crashworthiness in the automobile industry.

#### 2.2.6 Modified Mohr-Coulomb (MMC) fracture model

In the Mohr-Coulomb yield criterion, fracture is considered to occur when the shear stress (*τ*) and the vertical stress (*σ _{n}*) acting on the fracture surface reach a threshold value (

*c*

_{2}) (Eq. (24)). Bai and Weirzbicki (2010) proposed the Modified Mohr-Coulomb (MMC) model that expressed the existing Mohr-Coulomb fracture condition on the basis of strain.

*A*and

*n*are materials constants of the constitutive equation (Eq. (35)). This study did not treat the experiments on the Lode angle in negative domain. Therefore,

*c*

_{1},

*c*

_{2},

#### 2.2.7 Hosford-Coulomb (HC) model

After Mohr and Marcadet (2015) substituted the Tresca equivalent stress in the Mohr-Coulomb yield criterion with Hosford stress (Eq. (28)), they proposed the Hosford-Coulomb model for which the left hand side was changed from effective stress to effective plastic strain (see Eq. (30)). The HC model includes second principal stress (*σ*_{2}) and can consider the ratio of equi-biaxial stress (*σ*_{2}/*σ*_{1}) that the MMC model does not consider. The HC fracture model is composed of four material constants (*a, b, c, n _{f}*). Roth and Mohr (2016) proposed 0.1 as the value of

*n*about general steel. Therefore, there are only

_{f}*a, b, c*remaining as the material constants of the fracture model. The HC model has been used by many researchers in recent years because the number of model constants is smaller than that of other fracture models and the material constants can be determined with a minimum of experiments (Roth and Mohr, 2016; Erice et al., 2017; Park et al., 2019; Cerik et al., 2019b).

### 2.3 Damage evolution model

The damage evolution model expressed as the accumulation of equivalent plastic strain was used in order to consider loading path and stress path effects (Eq. (31)). Fig. 3 shows the basic development shapes of damage indicator (*m*) based on the damage indicator (*D*). When *m* = 1, the model can be expressed as Eq. (32) and is considered as linear damage development model. If *m* is bigger than 1, the damage indicator based on the increase of plastic strain is expressed in exponential form. If *m* is smaller than 1, the damage indicator is expressed in square root form. As materials and structures experience external forces, fracture is considered to occur when the damage indicator reaches 1.0. In this study, linear damage accumulation model (*m* = 1.0) is used.

## 3. Experiment and numerical analysis

### 3.1 Tensile experiment

The fracture experiment about EH36 steel conducted in the precedent research (Park et al., 2019) was used as the experimental model for the calibration about material constants of the fracture model. Therefore, a brief description of fracture experiment and numerical analysis performed in the precedent experiment is required. Fig. 4 shows the drawings and names of specimens of the precedent research. The base plate thickness of EH36 steel used for manufacturing specimens is 25 mm, and the specimens were manufactured at 2 mm from the middle layer of the base plate thickness direction. The chemical compositions are listed in Table 1. The flat bar specimens (FB) were manufactured in compliance with ASTM (2004) specifications. The notched tension specimens (NT) were manufactured with the radius (R) of 20 mm. The central hole specimens (CH) have a hole with the radius of 3 mm in their center and were designed to induce pure tensile status at the fracture point. The shear specimens (SH) were designed to trigger shear fracture. The fracture experiment was performed at 0.5 mm/min through the stroke displacement control and displacements and loads of 50 mm extensometer and load cell were measured.

### 3.2 Correction of flow stress material constant

The yield function (*f*) of isotropic materials is composed of von Mises weighing stress (*k*), which are developed from materials (Eq. (33)). The plastic potential theory uses the associated flow rule based on the flow rule. In other words, the plastic strain increment vector is parallel to the outward pointing unit normal vector perpendicular to a yield plane at the stress point.

Plastic hardening equation of materials is expressed as the Swift-Voce constitutive equation (Eq. (34)) that is a linear combination of Swift constitutive equation (Eq. (35)) and Vove constitutive equation (Eq. (36)). In addition, the flow stress was considered by considering the yield plateau that the initial yield stress is maintained and dividing the yield plateau into first half and latter half at *α*) was determined to enable the load-displacement curve obtained through numerical analysis about NT specimens to have experimental value and minim error. The determined material constants and flow stress curve are presented in Table 2 and Fig. 5 (Park et al., 2019).

### 3.3 Numerical analysis about tensile experiment

Fig. 6 shows finite element models of three specimens. The shear specimen performs symmetric modeling in the thickness direction. The notched specimen and the center hole specimen were subject to 1/8 modeling with the symmetry condition in the directions of thickness, length and width. In the numerical analysis, the tensile load was implemented by applying forced displacement to the upper node of each specimen. A commercial finite element analysis program, Abaqus/Explicit (Simulia, 2018) was used for numerical analysis, and 3-dimensional 8-node reduced integration element (C3D8R) was used. 10 notched elements of specimen were placed in the direction of thickness through convergence test.

It is impossible to measure the fracture starting point when any fracture occurs in a specimen even though the fracture experiment is observed by using optical equipment. Even if fracture occurs on the surface, it is difficult to find an exact fracture starting point since the propagation progresses rapidly after nucleation, growth and coalescence of void, which is considered as the generation of fracture. In this study, the fracture initiation displacement was regarded as the point where the load was reduced sharply in the experiment. In addition, the fracture occurrence point was determined to be the factor with the maximum equivalent plastic strain at fracture displacement. Fig. 7 shows the load-extensometer displacement of experiment and numerical analysis to the fracture point. It can be seen that the numerical analysis accurately simulates the experimental results.

Fig. 8 shows the loading path based on the increment of equivalent plastic strain. It can be seen that the loading path of the CH specimen keeps the uniaxial tension state (*η*=0.33, *η*=0.0,

## 4. Material constants of fracture models

To determine material constants of fracture models, the optimization technique was applied. In other words, the material constants of each fracture model were set as design parameters, and the range of the design parameter was used as constraints. In addition, the minimization of the residual sum of squares (*R*^{2}) was set as objective function. For optimization, the FMINUNC function of MATLAB was used. As shown in Fig. 8, the residual sum of squares indicates the sum of squares of ratio between the loading path data obtained from numerical analysis and the prediction data *i* is the number of experiments used to calibrate the material constants.

Table 3 presents constraints of design parameters, final material constants and residual sum of squares *R*^{2}. The closer the residual sum of squares using Eq. (37), the less error. It can be seen that the error rate of the C-L-O model is the largest, and B-W, MMC, L-H and HC models relatively recently proposed are estimated to show the best fracture prediction because they represents nearly the same residual sum of squares.

Fig. 9 shows fracture models in terms of the three dimensional fracture strain plane and the plane stress condition with the stress triaxiality and the load angle as functions. It can be seen that all of the seven fracture models are not significantly affected by the stress triaxiality in three dimensions. This can be confirmed from the J-C model, which has a sole function of the stress triaxiality, and it can be confirmed that the J-C model is close to the fracture strain of about 1.3. In the C-L-O model, the fracture strain increases sharply as it is close to the cut-off region. In addition, it can be seen that five fracture models (MSS, B-W, MMC, HC and L-H) are sensitive to the variability of the Lode angle parameter and that the sensitivity of the Lode angle parameter in MSS is greatest. The four fracture models (B-W, L-H, MMC, and HC) relatively recently proposed do not show significant differences in fracture strain levels as a whole. This can be confirmed more accurately under plane stress conditions. In addition, the fracture models except for the MSS model have similar strain values in the region of 1/3⩽*η*⩽2/3. Therefore, each fracture model is expected to have similar fracture prediction accuracy in tensile region (1/3⩽*η*⩽2/3). In the case of the MSS model, it can be seen that the fracture strains in uniaxial tension and equi-biaxial tension are more than twice as high as the other models.

Fig. 10 shows the damage by fracture model based on the development of the equivalent plastic strain at fracture points of each specimen along with the equivalent plastic strains at the fracture points. The damage development processes for three specimens in B-W, MMC, L-H and HC models are similar. It can be seen that the C-L model has high accuracy about fracture prediction in the fracture experiments (CH and NT) on the tensile region, but damage indicator is under-estimated at the fracture point. The MSS model has significantly low fracture prediction accuracy about the CH specimen, which maintains the loading path in the uniaxial tension state.

## 5. Conclusion

In this paper, we presented the fracture behavior characteristics about steel grade EH 36 through the examination of various existing fracture models. Various phenomenological fracture models were represented as stress triaxiality and Lode angle parameter, which were known as dominant variables in fracture behavior of ductile materials. In addition, the flow stress of EH36 steel was extracted by deriving the material constants of the Swift-Voce constitutive equation through the numerical analysis about the experiments performed in the precedent research (Park et al., 2019). It could be seen that the flow stress, which had presented through numerical analysis about each experiment, exactly analogized the behavior of the steel to the large strain section. We drew the loading path at the fracture point and corrected the material constants of each fracture model. It could be seen that the fracture model was illustrated in stress triaxiality, Lode angle and fracture strain spaces and that EH36 steel was not sensitive to the effect of the stress triaxiality. In addition, It was confirmed that the B-W, L-H, NNC, HC models, which has been relatively recently proposed fracture models, had similar fracture strain levels.

The fracture models examined in this study is based on solid element. For the fracture simulation of structures with thin plate compared to their size (for example, ship structures and marine structures), the use of shell element is essential due to time constraints for modeling and analysis. Recently, multiple researchers has proposed various methods to derive results approximate to solid element in shell element based fracture simulation (Walters et al., 2014; Kõrgesaar et al., 2014; Park and Mohr, 2017). To apply the fracture models to the actual marine structures, more studies on shell element based fracture simulation technique are required.

Steel grade used in industry is classified according to operating guarantee temperature, minimum yield strength, impact toughness and others. The minimum yield strength and tensile strength of high tensile strength steel with Grade AH36, DH36 and EH36 specified in DNVGL regulation (DNVGL, 2018) are 355 MPa and 490 MPa, respectively, and the upper limits of yield strength and fracture elongation ratio are not clearly specified. Therefore, steels with the same grade may show the differences in mechanical properties. The material hardening behavior directly affects the loading path, and reducing the uncertainty of the material hardening behavior is an important issue for increasing accuracy in predicting the fracture behavior. Cho et al. (2015) statistically presented the yield strength and the tensile strength of various steel grades (Table 4). The yield strength and tensile strength of the EH36 steels examined in this study are 428.028 MPa and 523.2 MPa, respectively, which are not significantly different from the average value of EH36 steel given in Table 4. Therefore, it is considered that the fracture model of EH36 steel can represent the fracture model of the same grade steel. However, it can be seen that AH36 and DH36 steels have higher tensile strength than EH36 steel. The construction of database through the fracture experiments on various steels is essential in order to make general conclusions about fracture behavior of steel.

## Postscript

This study is the research result performed under the support of Korea Research Institute of Ships & Ocean Engineering's ‘Study on the basis of structural engineering design for the construction of structural safety evaluation system (3/3) [PES3250]’, National Research Foundation of Korea’s ‘Overseas excellent researcher support project (2017H1D3A1A01055137)’ and Ministry of Trade, Industry and Energy’s ‘Korea and UK offshore plant global expert training project’.