Anisotropic Fracture Limit for the Failure Prediction of ...

1 downloads 0 Views 21MB Size Report
Mar 31, 2017 - Continuum theory of plasticity. John Wiley & Sons, 1995. ..... [194] Y. Zhu, B. Dodd, R. Caddell, and W. Hosford. Convexity restrictions on non- ...
Anisotropic Fracture Limit for the Failure Prediction of Advanced High-Strength Steel Sheets

Park, Namsu Advisor: Huh, Hoon

A dissertation submitted to the faculty of Korea Advanced Institute of Science and Technology in partial fulfillment of the requirements for the degree of Doctor of Philosophy in Mechanical Engineering

Daejeon, Korea March 31, 2017 Approved by

Huh, Hoon Professor of Mechanical Engineering

The study was conducted in accordance with Code of Research Ethics1 . 1

Declaration of Ethical Conduct in Research: I, as a graduate student of Korea Advanced Institute of Science

and Technology, hereby declare that I have not committed any act that may damage the credibility of my research. This includes, but is not limited to, falsification, thesis written by someone else, distortion of research findings, and plagiarism. I confirm that my thesis contains honest conclusions based on my own careful research under the guidance of my thesis advisor.

DME 20135111

박남수. 고강도 강판의 파괴 예측을 위한 이방성 파괴 한계에 관한 연 구. 기계공학과. 2017년. 182+xiv 쪽. 지도교수: 허훈. (영문 논문) Park, Namsu. Anisotropic Fracture Limit for the Failure Prediction of Advanced High-Strength Steel Sheets. Department of Mechanical Engineering. 2017. 182+xiv pages. Advisor: Huh, Hoon. (Text in English)

초록 최근 자동차 산업 분야에서 차체의 충돌성능과 안전성을 위하여 차체에 사용되고 있는 재료 중 고강도강의 비율을 점차적으로 확대하고 있다. 그러나 고강도강(Advanced High Strength Steels, AHSSs)은 기계적 강도가 높은 반면 연성이 작기 때문에 성형공정에서 적은 량의 소성 변형이 수반되는 경우에도 재료의 파단이 발생하는 현상이 관찰되었다. 네킹을 거의 동반하지 않고 파단이 발생하므로 기존의 네킹 기반의 성형한계도(Forming Limit Diagram, FLD)는 파 단을 예측하기에는 한계가 있다. 파괴가 발생하는 하중경로 역시 기존 성형한계도가 다루지 않는 단축 압축에서 단축인장 사이의 하중경로에서 파괴가 발생하기 때문에 다양한 하중경로에 대한 파괴기반의 성형성 평가가 요구되고 있는 실정이다. 판재는 압연공정으로 인하여 일반적으로 두께방향 및 평면 상의 이방성을 가지며 이는 하중방향 및 하중경로와 더불어 재료의 파괴에 영향을 미친다. 고강도강을 이용한 차체 성형 시 재료 성형성의 정량적 평가를 위해서는 판재의 이방성이 파괴물성에 미치는 영향을 필수적 으로 고려하여야 한다. 따라서, 본 연구에서는 다양한 하중경로 및 하중방향을 고려한 연성파 괴조건을 개발하고 이를 바탕으로 자동차 산업의 파괴조건 적용을 위한 여러 파괴성형한계를 제시하고자 한다. 기존의 Lou–Huh 연성파괴조건을 기반으로 Hill’s 48 이방성 항복조건을 적용 하여 판재의 파괴에 영향을 미치는 이방성을 고려할 수 있도록 하였고, 코시 응력 텐서(Cauchy stress tensor)의 좌표 변환을 통해 주응력 및 최대주응력 방향에 의존적인 이방성 파괴조건을 제안하였다. 제안된 파괴조건을 정량적으로 평가하기 위해 고강도강의 일종인 DP980 1.2t 판재 내 다양한 하중방향(rollding direction, diagonal direction, transver direction)에서 획득한 파괴물성을 이용하여 파괴조건식의 계수를 결정하였고, 개발된 파괴조건을 통하여 각 방향 별 파괴물성을 정량화하였다. 제안된 파괴조건은 로데 변수(Lode parameter)와 Hill’s 48 이방성 항복조건 기반의 응력삼 축성 및 최대 주응력 방향으로 정의되며, 이 파라미터들의 조합을 통해 재료점이 겪는 다양한 하중경로를 고려할 수 있다. 또한, 판재의 평면 내 이방성이 특정 하중경로와 더불어 파단물성 에 미치는 복합적인 효과가 최대 주응력 방향으로 연관되어 판재의 평면 내 이방성을 고려한 특정하중경로에서의 파단변형률을 정량적으로 평가할 수 있다. 사각 컵 드로잉 시험과 특정 하중경로를 구현한 시편의 인장시험 결과 및 제안된 파괴조건의 파괴예측결과의 비교를 통해 파괴조건의 파괴예측성능을 확인하였고, 고강도강의 파괴예측을 위한 이방성 연성파괴조건으

i

로써의 유효성을 검증하였다. 핵 심 낱 말 이방성, 변형률 기반 파괴성형한계도, 응력 기반 파괴성형한계도, 극좌표계 기반 파괴성형한계도, 이방성 연성파괴조건

Abstract This thesis is concerned with developing anisotropic fracture forming limit criteria for the prediction of the material formability in sheet metal forming to predict the sudden fracture in complicated forming processes for advanced high-strength steel (AHSS) sheets. A new anisotropic ductile fracture criterion is developed to provide a phenomenological model in consideration of non-directionality of the equi-biaxial fracture strain as well as the effect of anisotropy at the macroscopic level particularly for sheet metal forming application. In developing the anisotropic ductile fracture criterion, the Hill’s 48 criterion is employed to take account of the influence of anisotropy on the equivalent plastic strain at the onset of fracture. For the derivation of the anisotropic ductile fracture criterion, the principal stresses (σ1 , σ2 , σ3 ) are expressed in terms of the stress triaxiality, the Lode parameter, and the effective stress (ηH , LP , σ ¯H ) based on the Hill’s 48 criterion to include the effect of directionality of typical loading states on the material orientation. For determination of the parameters of the fracture criterion proposed, the two-dimensional Digital Image Correlation (DIC) method is utilized to measure the strain histories on the surface of three different types of specimens: the pure shear; the uniaxial tension; and the plane strain tension. The measurement results are investigated to quantitatively identify the anisotropy effect on the equivalent plastic strain at the onset of fracture. From the proposed fracture criterion, three different kinds of fracture forming limit criteria are introduced: a strain-based Fracture Forming Limit Diagram (FFLD); a stress-based FFLD; and a Polar Effective Plastic Strain (PEPS) FFLD. These fracture forming limit criteria are obtained with an assumption of the proportional loading under the plane stress condition. A scaling method for a strain-based fracture forming limit criterion is also discussed in order to capture the onset of fracture using a single forming limit curve for an anisotropic material. Experimental validations are performed in the viewpoint of a structural and a specimen level through the square cup drawing and the tensile test at additional loading directions. In comparison of the experimental results with the ones predicted from the proposed fracture criterion, it is clearly concluded that the proposed fracture criterion has a sufficient capability to predict the fracture initiation over a wide range of stress states in consideration of the material anisotropy as well as non-directionality of the equi-biaxial fracture strain. Comparison of the

ii

strain-based FFLD and the PEPS FFLD reveals that the PEPS FFLD is able to predict the fracture initiation regardless of the strain path change whereas the performance of the strainbased FFLD can be valid only for the case that the material point is under the proportional loading condition. Consequently, the PEPS FFLD can play a key role in predicting the fracture initiation during complicated sheet metal forming processes of AHSS sheets in addition to its path independence and simplicity of measuring strains in real forming processes. Keywords Anisotropy, Strain-based Fracture Forming Limit Diagram (Strain-based FFLD), Stress-based Fracture Forming Limit Diagram (Stress-based FFLD), Polar Effective Plastic Strain Fracture Forming Limit Diagram (PEPS FFLD), Anisotropic Ductile Fracture Criterion.

iii

Contents Abstract

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

i

Contents

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

iv

List of Tables

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

vii

List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

ix

Chapter 1

Introduction

1

1.1

Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1

1.2

Literature review . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

1.2.1

Ductile fracture criteria . . . . . . . . . . . . . . . . . . .

2

1.2.2

Modeling of anisotropic yield behavior . . . . . . . . . .

5

Objective and scope of the thesis . . . . . . . . . . . . . . . . . .

8

1.3

Chapter 2 2.1 2.2

Theoretical Backgrounds

9

Characterization of the material anisotropy . . . . . . . . . . . .

9

2.1.1

Hill’s 48 criterion . . . . . . . . . . . . . . . . . . . . . . .

9

Anisotropic coefficients calculation . . . . . . . . . . . . . . . . .

12

2.2.1

Calculation of the yield stress in an arbitrary direction

13

2.2.2

Calculation of the r-value in an arbitrary direction . . .

13

Chapter 3

Quantification of Plastic Anisotropy at Various Orientations

16

3.1

Tensile tests at various orientations . . . . . . . . . . . . . . . . .

16

3.2

Mesurement of the r-value at various orientations . . . . . . . .

17

Chapter 4 4.1

Characterization of the Stress State

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηv , θ, σ ¯v ) . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2

24

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηH , LP , σ ¯H )

4.4

21

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηv , LP , σ ¯v ) . . . . . . . . . . . . . . . . . . . . . . . . . .

4.3

21

. . . . . . . . . . . . . . . . . . . . . . . . .

27

Lode dependence in the principal strain space for anisotropic materials . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

iv

30

Chapter 5 5.1

Development of an Anisotropic Ductile Fracture Criterion

Construction of an uncoupled anisotropic ductile fracture criterion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.2

33 33

Determination of the parameters of the proposed fracture criterion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Chapter 6

Construction of Fracture Forming Limit Diagrams

6.1

Strain-based fracture forming limit diagram . . . . . . . . . . . .

6.2

Scaling of the strain-based fracture forming limit diagram ac-

36 41 41

cording to the direction of major stress axis . . . . . . . . . . . .

42

6.3

Stress-based fracture forming limit diagram . . . . . . . . . . . .

44

6.4

Polar effective plastic strain-based fracture forming limit diagram 46

Chapter 7

Emprical Way of Modeling the Fracture Criterion Considering the Material Anisotropy

7.1

48

Model parameters dependent on the maximum principal stress direction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.2

Weight function considering non-directionality of the equi-biaxial fracture strain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Chapter 8

8.2

sponse in the large deformation range . . . . . . . . . . . . . . .

52

8.1.1

Theory for stress and strain measurement . . . . . . . .

52

8.1.2

Flow stress curves at the three loading directions considering large deformation . . . . . . . . . . . . . . . . . .

55

Fracture limits under the proportional loading condition . . . .

56

Prediction of fracture strains at additional loading directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Chapter 9

9.2

51

Hydraulic bulge tests for quantification of the stress–strain re-

8.2.1

9.1

49

Discussion of Anisotropic Fracture Forming Limits for Application to the Automotive Industry

8.1

48

Numerical Implementation

57 59

Flow rule . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

59

9.1.1

Non-associative flow rule . . . . . . . . . . . . . . . . . . .

61

Explicit time integration . . . . . . . . . . . . . . . . . . . . . . .

62

9.2.1

64

Chapter 10

Stress update algorithm . . . . . . . . . . . . . . . . . . . Experimental Validation

10.1 Square cup drawing test

67

. . . . . . . . . . . . . . . . . . . . . . . v

67

10.2 Square cup drawing simulation

. . . . . . . . . . . . . . . . . . .

67

10.3 Comparison of the experimental and simulation results . . . . .

69

Chapter 11

Conclusions and Future Works

71

11.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

71

11.2 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

72

Appendix A.

Generalization of the Anisotropic Stress Triaxiality

73

A.1 Anisotropic stress triaxiality defined by the Yld2004-18p criterion 73 A.2 Construction of the uncoupled anisotropic fracture criterion based on the Yld91 criterion . . . . . . . . . . . . . . . . . . . . . Appendix B.

78

Construction of the Anisotropic Fracture Forming Limit Diagram Based on the Non-quadratic Yield Criterion

Bibliography

84 87

Acknowledgments in Korean

183

Curriculum Vitae in Korean

185

vi

List of Tables Table 1

Uncoupled isotropic ductile fracture criteria . . . . . . . . . . . . . . . . . . 104

Table 2

Uncoupled anisotropic ductile fracture criteria . . . . . . . . . . . . . . . . . 105

Table 3

Material quantities needed for the identification of several yield criteria [12] 106

Table 4

Chemical compositiona of the DP980 1.2t steel sheet, in wt% . . . . . . . . 107

Table 5

Mechanical properties of the DP980 1.2t at the three loading directions . . . 107

Table 6

Parameters of the mixed Swift–Voce model of the DP980 1.2t steel sheet at

the three loading directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 Table 7

Equivalent plastic strains to fracture evaluated at the three loading directions107

Table 8

Typical anisotropy coefficients for various materials [76]

Table 9

Parameters of the asymptotic regression model at the three loading directions108

. . . . . . . . . . . 108

Table 10 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the r-values obtained in an average sense

. . . . . . . . . . . . . . . . . . . . . . 108

Table 11 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the instantaneous r-values after saturated . . . . . . . . . . . . . . . . . . . . . . 108 ¯ ξ, χ and LP Table 12 Typical stress condition with corresponding values of θ, θ, (σ1 ≥ σ2 ≥ σ3 ) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 Table 13 Combination of the stress triaxiality and the Lode parameter denoting typical stress states . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 Table 14 Combination of the anisotropic stress triaxiality and the Lode parameter for the DP980 1.2t denoting the typical stress states where the parameters of Hill’s 48 criterion are obtained by the r-value evaluated in an average sense . . . . . . 110 Table 15 Model parameters of the anisotropic ductile fracture criterion . . . . . . . . 110 Table 16 Model parameters of the anisotropic ductile fracture criterion constructed in a purely empirical way . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110 Table 17 Anisotropic characteristic parameters of the Hill’s 48 criterion determined by the normalized yield stresses where the equi-biaxial yield stress is set as the yield stress for the rolling direction . . . . . . . . . . . . . . . . . . . . . . . . . . 110 Table 18 Conditions for the hydraulic bulge test of the DP980 1.2t . . . . . . . . . . . 110

vii

Table 19 Parameters of the mixed Swift–Voce model of the DP980 1.2t steel sheet in the large deformation range . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 Table 20 Mechanical properties of the DP980 1.2t steel sheet evaluated from the hydraulic bulge test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 Table 21 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the saturated values of the normalized yield stress ratios . . . . . . . . . . . . . . 111 Table 22 Model parameters of the anisotropic ductile fracture criterion constructed under the proportional loading condition . . . . . . . . . . . . . . . . . . . . . . . 111 Table 23 Equivalent plastic strains to fracture evaluated at the additional loading directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 Table 24 Numerical algorithm for the tangent cutting plane integration method . . . 112 Table 25 Conditions for the square cup drawing test of the DP980 1.2t steel sheet . . 113 Table 26 Specifications of 3,000 kN hydraulic servo press . . . . . . . . . . . . . . . . 113 Table 27 Total number of elements and nodes generated for each FE model . . . . . . 113 Table 28 Parameters of the exponential fitting function for each simulation group . . 113 Table 29 Anisotropic characteristic parameters of the Yld91 criterion determined by the normalized yield stress ratios saturated . . . . . . . . . . . . . . . . . . . . . 113

viii

List of Figures Figure 1

Schematic forming limits classified by: (a) Necking, wrinkling and shear

fracture (after Martins et al. [127]); (b) Deformation modes of the circular grid . 114 Figure 2

Fracture of Dual-Phase (DP) steel sheets in various sheet metal forming

processes (after research works in ductile fracture of advanced high-strength steel sheets [1, 112, 124, 155, 179, 190]) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 Figure 3 Nucleation, growth, and coalescence of voids during the tension test (after ˇ Sebek [154]) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 Figure 4 Drawings of test specimens: (a) Pure shear specimen; (b) Dog-bone specimen; (c) Plane strain grooved specimen [mm] . . . . . . . . . . . . . . . . . . . . 117 Figure 5 Uniaxial tensile testing machine (INSTRON5583) . . . . . . . . . . . . . . . 118 Figure 6 True stress–true strain curves of the DP980 1.2t steel sheet evaluated at the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) Flow stress curves plotted all together . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 Figure 7 Strain paths evaluated at the three loading directions (0°, 45°, 90°) under typical loading conditions: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 Figure 8

Equivalent plastic strains at the onset of fracture evaluated at the three

loading directions (0°, 45°, 90°) under typical loading conditions: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case; (d) Comparison among the ones averaged

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120

Figure 9 Representative strain paths evaluated at the three loading directions (0°, 45°, 90°) with typical loading conditions . . . . . . . . . . . . . . . . . . . . . . . 121 Figure 10 Transverse strain versus longitudinal strain distribution of specimens at the three loading directions: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . 122 Figure 11 Instantaneous r-values fitted by the asymptotic regression model and the corresponding residuals: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . 123 Figure 12 Transverse plastic strain versus longitudinal plastic strain distribution predicted by the intantaneous r-values fitted by the asymptotic regression model and the corresponding residuals: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . 124

ix

Figure 13 Instantaneous normal and plannar anisotropic quantities according to the plastic work and the equivalent plastic strain . . . . . . . . . . . . . . . . . . . . 125 Figure 14 Comparison of the normalized yield loci constructed with Hill’s 48 and von Mises criteria . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 Figure 15 Representation of the principal stress state in Haigh–Westergaard space (σ1 , σ2 , σ3 ) and on the deviatoric plane (s1 , s2 , s3 ) . . . . . . . . . . . . . . . . . . 126 Figure 16 Dependence of: (a) The stress triaxiality on various indicators concerning the deviatoric stress state; (b) Various state variables on the normalized third invariant of the deviatoric stress tensor . . . . . . . . . . . . . . . . . . . . . . . . 127 Figure 17 Rotation axis and angles for the coordinate transformation . . . . . . . . . . 128 Figure 18 Stress states in the space of the newly defined stress triaxiality ηH and the Lode parameter LP according to the three loading directions of 0°, 45°, and 90°: (a) The change of new combination of the ηH and LP for DP980 1.2t steel sheet; (b) The change of new combination of the ηH and LP for SPCC 1.0t steel sheet; (c) The ηH corresponding to the typical stress paths for SPCC 1.0t steel sheet; (d) The LP corresponding to the typical stress paths . . . . . . . . . . . . . . . . 130 Figure 19 Cut-off region: (a) The initial one; (b) The modified one . . . . . . . . . . . 131 Figure 20 Strain-based 3D fracture envelopes according to the three loading directions: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 Figure 21 Fracture loci according to the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) Fracture loci plotted all together . . . . . . . . . . . . . . . . . . . . . . . 133 Figure 22 Damage evolutions evaluated according to the loading paths and the three loading directions: (a) RD; (b) DD; (c) TD . . . . . . . . . . . . . . . . . . . . . 134 Figure 23 Fracture forming limit diagram (FFLD) according to the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) FFLDs plotted all together . . . . . . . . . . 135 Figure 24 Fracture loci and its corresponding FFLDs at different directions of the maximum principal stress: (a) Fracture loci plotted according to the maximum principal stress direction; (b) FFLDs plotted according to the maximum principal stress direction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136 Figure 25 Major fracture strain in plane strain according to the maximum principal stress direction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137 Figure 26 Scaling factor in the space of the maximum principal stress direction and the strain path . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137

x

Figure 27 Fracture forming limit diagrams at different directions of the maximum principal stress: (a) FFLDs scaled at different directions of the maximum principal stress; (b) Comparison between FFLDs predicted from the anisotropic ductile fracture criterion and the FFLDs scaled . . . . . . . . . . . . . . . . . . . . . . . 138 Figure 28 Stress-based 3D fracture envelopes according to the three loading directions: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 Figure 29 Stress-based fracture forming limit criteria according to the three loading directions (0°, 45°, 90°): (a) The stress-based FFLDs; (b) The corresponding fracture loci . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140 Figure 30 Stress-based fracture forming limit criteria according to various directions of the maximum principal stress: (a) The stress-based FFLDs; (b) The corresponding fracture loci . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141 Figure 31 Stress-based fracture forming limit diagram: (a) Not considering the one point convergence at the equi-biaxial fracture stress; (b) Considering the one point convergence at the equi-biaxial fracture stress . . . . . . . . . . . . . . . . . 142 Figure 32 Polar Effective Plastic Strain (PEPS) Fracture Forming Limit Diagrams (FFLDs) according to the three loading directions (0°, 45°, 90°) . . . . . . . . . . 143 Figure 33 Comparison of the Polar Effective Plastic Strain (PEPS) Fracture Forming Limit Diagrams (FFLDs) predicted from the fracture criterion with the FFLDs scaled from the strain-based FFLD for RD . . . . . . . . . . . . . . . . . . . . . . 144 Figure 34 Damage evolutions according to the loading paths and three loading directions when the empirical way is employed to the newly developed fracture criterion: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . . . . . . . . . . 145 Figure 35 Fracture forming limit diagram (FFLD) according to the three loading directions when the empirical way is employed to the newly developed fracture criterion: (a) 0°; (b) 45°; (c) 90°; (d) FFLDs plotted all together

. . . . . . . . . 146

Figure 36 Schematic description of the radius of the curvature . . . . . . . . . . . . . . 147 Figure 37 Hydraulic bulge test system in POSCO [7] . . . . . . . . . . . . . . . . . . . 147 Figure 38 Flow stress curves obtained from the hydraulic bulge tests . . . . . . . . . . 148 Figure 39 Hydraulic bulge test result: (a) Fractured specimen; (b) Equivalent plastic strain distribution before the onset of fracture . . . . . . . . . . . . . . . . . . . . 149 Figure 40 Schematic description of the plastic work equivalence principle . . . . . . . . 150

xi

Figure 41 Plastic works per volume of the DP980 1.2t steel sheet and the directional yield stresses normalized by the equi-biaxial yield stress: (a) Plastic works per volume; (b) Normalized yield stress ratios . . . . . . . . . . . . . . . . . . . . . . 151 Figure 42 Flow stress curves in the large deformation range . . . . . . . . . . . . . . . 152 Figure 43 Strain-based 3D fracture envelopes constructed under the proportional loading condition: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . . . . . . . . . . . . . . . . 153 Figure 44 Fracture loci constructed under the proportional loading condition: (a) Predicted from the fracture criterion considering the model parameters dependent on the loading direction; (b) Predicted from the fracture criterion considering the model parameters dependent on the loading direction in combination with the weight function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 154 Figure 45 Strain-based FFLDs constructed under the proportional loading condition predicted from the fracture criterion considering the model parameters dependent on: (a) The loading direction; (b) The loading direction in combination with the weight function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 155 Figure 46 Comparison between anisotropic and isotropic fracture forming limits: (a) strain-based FFLDs with the necking-based FLD; (b) fracture loci . . . . . . . . 156 Figure 47 Fracture forming limit envelopes constructed under the proportional loading condition: (a) Strain-based fracture forming envelope; (b) Stress-based fracture forming envelope . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157 Figure 48 Three-dimensional fracture envelopes constructed under the proportional loading condition: (a) In the principal strain space; (b) In the principal stress space . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158 Figure 49 Comparison of the fracture strains obtained from the experiments with the ones predicted from the proposed fracture criterion under the proportional loading condition: (a) Pure shear; (b) Uniaxial tension; (c) Plane strain tension . . . 159 Figure 50 Prediction of the fracture strain over a wide range of stress state under the plane stress condition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160 Figure 51 Normal directions of the yield surface and the plastic potential in the nonAFR approach . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160 Figure 52 Relationship between the acceleration, velocity, and the displacement in the explicit time integration scheme . . . . . . . . . . . . . . . . . . . . . . . . . . . . 161 Figure 53 Graphical description of the tangent cutting plane algorithm . . . . . . . . . 161 Figure 54 Specimen fractured during the square cup drawing test . . . . . . . . . . . . 162 xii

Figure 55 Schematic descriptions of each square cup drawing simulation: (a) Symmetric case; (b) Non-symmetric case with a bias of 5 mm

. . . . . . . . . . . . . . . 163

Figure 56 Yield loci and plastic potential of the DP980 1.2t steel sheet . . . . . . . . . 164 Figure 57 Flow stress curves used in the FE analysis . . . . . . . . . . . . . . . . . . . 164 Figure 58 Finite elemement modeling of each specimen: (a) Pure shear; (b) Uniaxial tension; (c) Plane strain tension . . . . . . . . . . . . . . . . . . . . . . . . . . . . 165 Figure 59 Comparison of the load–displacement curves for each specimen: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case . . . . . . . . 166 Figure 60 Results of the FE analysis for a quarter FE model according to various mesh size . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 167 Figure 61 Results of the FE analysis for a half FE model according to various mesh size167 Figure 62 Crack predicted in FE simulation: (a) A quarter FE model; (b) A half FE model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 168 Figure 63 Localized equivalent plastic strain distribution before the fracture initiation with different mesh sizes of: (a) 0.5 mm; (b) 3 mm . . . . . . . . . . . . . . . . . 169 Figure 64 Punch stroke to fracture for each FE group: (a) A quarter FE model; (b) A half FE model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 170 Figure 65 Fracture prediction using the strain-based FFLD for an: (a) Outer layer; (b) Inner layer . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171 Figure 66 Fracture prediction using the strain-based FFLD according to the loading direction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172 Figure 67 Fracture prediction of an outer layer using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence 173 Figure 68 Fracture prediction of an inner layer using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence 174 Figure 69 Fracture prediction according to the loading direction using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 175 Figure 70 Fracture prediction using the PEPS FFLD at the punch edge of the sheet . 176 Figure 71 Fracture prediction using the PEPS FFLD at the punch radius of the sheet 176 Figure 72 Fracture prediction using the PEPS FFLD at the punch edge of the sheet where no fracture initiates . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 177

xiii

Figure 73 Conditions on the material parameters to reduce: (a) The Yld2004-18p criterion to the Yld91 criterion; (b) The Yld91 criterion to Hosford, Hill’s 48, von Mises, and Tresca criteria . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 178 Figure 74 Yield loci of AA6111-T4 constructed by the generalized anisotropic stress triaxiality under the plane stress condition . . . . . . . . . . . . . . . . . . . . . . 179 Figure 75 Rotation axis and angles defined for the application of the fracture surface to sheet metal forming . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 179 Figure 76 Distribution of the scaling values over a wide range of the stress states according to the three loading directions: (a) 0°; (b) 45°; (c) 90° . . . . . . . . . . 180 Figure 77 Scaling values according to the Lode parameter at the three loading directions181 Figure 78 Influence of the exponent m of the Yld91 criterion on the yield locus . . . . 181 Figure 79 Influence of the exponent m of the Yld91 criterion on the fracture forming limit: (a) Fracture loci; (b) Fracture forming limit diagrams . . . . . . . . . . . . 182

xiv

Chapter 1 Introduction 1.1

Background In the automotive industries, it is important to evaluate the forming severity of sheet metals

to make an efficient design of structural components for the auto-body. In dealing with the forming severity, necking has been regarded as a major failure mechanism because it gives the direct information on the loss of a load carrying capability. Various analytical models based on necking, thus, have been proposed so as to evaluate the forming limit of sheet metals accurately: the Hill’s localized necking model [87]; the Swift’s diffuse necking model [169]; the imperfectionbased Marciniak–Kuczynski model [125]; the modified maximum force criterion [89]; and the vertex theory [160,193]. Meanwhile, the experimental method to measure the forming limit was firmly established in the 1960s by Keeler and Backofen [96] and Goodwin [74]. This experimental methodology has been widely accepted as a standard tool and named as Forming Limit Diagram (FLD) for indicating the forming limit of sheet metals. The conventional FLD describes necking limits and thickness reductions in terms of the major and minor strains over the deformation modes from the uniaxial tension to the equi-biaxial tension as shown in Figure 1. Recently, the usage of advanced high-strength steels (AHSS) for automotive applications has been remarkably increasing due to a fact that AHSS sheets, in general, have outstanding strength to weight ratio, which allows lightweight construction of the car body structures with the use of thinner gages and improves the performance of crash relevant parts; however, AHSS sheets usually comes at expense of ductility resulting in having relatively poor formability compared with that of the conventional steel sheets [80]. With increasing demands of advanced high strength steels (AHSSs) including aluminum alloys and magnesium alloys to the automotive industry, there arises a challenging issue to deal with a sudden fracture of sheet metals during forming processes as represented in Figure 2. This undesired failure mode, which is often called the shear fracture, mainly results from the low ductility of advanced metal sheets and it shows the fracture surface slanted along the maximum shear stress direction through the thickness of metal sheets with the small amount of necking. In addition, this failure phenomenon is observed not only in tension but also in shear and compression conditions where the thickness reduction is negligible. Underlying this challenging issue is a limitation to accurately predict the onset of fracture in aid of the necking-based conventional FLDs in a vehicle development

1

stage because it is not able to deal with the sudden fracture in high strength steel sheets, which lies in a fact that AHSS sheets generally fail with little amount of necking different from that with the conventional steel sheets. Meanwhile, ductile fracture in advanced metals is observed in shear and compression conditions where the value of the stress triaxiality is low or negative in bulk metal forming as reported by Bao and Wierzbicki [17], Børvik et al. [29] and Khan and Liu [98,99]. From these experimental observations, ductile fracture criterion can be regarded as a potential solution to deal with the sudden fracture of AHSS sheets instead of a necking-based forming limit diagram which neither monitors the region for low and negative stress triaxiality nor predicts the sudden fracture; in-depth analysis and modeling of fracture phenomena are essential for accurate predictions with the aid of Computer-Aided Engineering (CAE). In general, sheet metals exhibit a significant anisotropy of the mechanical properties due to their crystallographic structure and the characteristics of the rolling process [11]. This change of crystallographic structure in the material has a notable influence on the fracture behavior [139]. Consequently, there remains a strong need for experimental examination and the modeling of ductile fracture criterion to effectively and efficiently deal with the challenging issue considering the material anisotropy for the application of AHSS sheets to the automotive industry especially for sheet metal forming processes.

1.2

Literature review

1.2.1

Ductile fracture criteria

The failure of advanced metal sheets can be considered as the ductile fracture with some amount of plastic deformation since the onset of fracture comes after the concurrent and mutual interactive fracture mechanisms with regard to the nucleation, growth, and coalescence of microcavities or voids. The forming severity of advanced metal sheets is, therefore, possible to be assessed using the ductile fracture criterion instead of the necking-based forming limit criterion. The reason behind this lies in a fact that the ductile fracture can take place over a wide range of stress states including shear and compression conditions under low and negative stress triaxiality. As illustrated in Figure 3, the ductile fracture usually follows a multi-step failure process from the micro-mechanical point of view such as: nucleation of the microscopic voids by either fracture or decohesion of the second-phase particles and inclusions; growth of the microscopic voids due to the localized plastic deformation; localization of the plastic flow between the growing voids; final tearing or rupture of the ligaments between the grown and enlarged voids [159,172]. Over the past decades, numerous approaches have been proposed in the literature that

2

relies on various features to deal with the onset of fracture at a wide range of stress states. The proposed approaches can be classified into two categories such as “coupled” and “uncoupled” fracture criteria. The main difference between these two categories of fracture criteria is whether the strength of metals is affected by damage accumulated through the nucleation, the growth, and the coalescence of voids. Coupled fracture criteria are firmly established based on the postulate explained above; on the contrary, uncoupled fracture criteria do not consider the effect of accumulated damage on the load carrying capability of metals. In the case of coupled fracture criteria, the nucleation and the growth of voids are mathematically formulated based on porous plasticity introduced by Gurson [78], which lies in a fact that the stress triaxiality is the main parameter that controls the void growth [128, 150]. This approach has gained great attention after Tvergaard and Needleman [174] incorporated the effect of void coalescence, which is named as the Gurson–Tvergaard–Needleman (GTN) ductile fracture criterion. In the GTN ductile fracture criterion, the accumulated damage is represented by the void volume fraction which is coupled by the constitutive equation so as to induce the effect of strength weakening. As a base formula, the GTN ductile fracture criterion has been widely utilized to predict shear localization at the low stress triaxiality [133, 135, 184] as well as void shape changes [38, 56, 73]. Another popular coupled fracture criteria are based on the Continuum Damage Mechanics (CDM) introduced by Kachanov [94], which is further improved by various research works [34–37, 43–45, 109–111, 153]. The CDM considers the damage as an internal variable which grows with plastic deformation to represent the local distribution of micro-defects which eventually induce the change of material properties such as macroscopic elastic and hardening modulus. In the case of uncoupled fracture criteria, dozens of approaches have been proposed based on microscopic mechanisms and experimental observations of the ductile fracture with various hypotheses: Freudenthal [68]; Cockcroft and Latham [50]; Brozzo et al. [33]; Oh et al. [136]; Oyane et al. [138]; Clift et al. [49]; and Ko et al. [103] and so on. These criteria were widely applied to deal with various engineering problems such as compressive upsetting tests, axisymmetric extrusion, strip compression and tension due to their simple formulae with a few parameters for experimental calibration. Uncoupled fracture criteria reviewed above are, however, incapable of fully describing the ductile fracture behavior over a wide range of stress states due to their hypotheses for solving special problems. Recently, Bao and Wierzbicki [17] proposed a phenomenological fracture criterion based on results of extensive fracture tests of AA2024-T351 including upsetting, shear, and tensile tests, where a wide range of stress state was covered. Their experimental results revealed an important clue that the stress triaxiality is not the only 3

factor that affects the ductile fracture by showing the non-monotonic dependence of the stress triaxiality on the fracture strain. Studies by Kim et al. [100, 101] and Gao et al. [69] on the void-contained cell model found that the Lode parameter should be taken into consideration so as to distinguish the stress state with the same stress triaxiality because the macroscopic stress–strain response, as well as the void growth and coalescence behavior of a representative material volume significantly, varies for each stress state even though the stress triaxiality stays the same. Lode angle dependence on ductile fracture was reported by Xue [182] and Xue and Wierzbicki [185] and it was confirmed that the Lode angle parameter, which is related to the third invariant of the stress deviator tensor, and the stress triaxiality have a mutual influence on the initiation of ductile fracture from various experimental results performed by Barsoum and Faleskog [24, 25] and Korkolis and Kyriakides [105]. In their experiments, a shift in rupture mechanism from void internal necking to internal shearing takes place as the stress triaxiality decreases with an increase in the amount of shearing, which is accompanied with not only a dramatic change in the Lode parameter but also a significant drop in ductility. Bai and Wierzbicki [10] modified the classical Mohr–Coulomb criterion which was basically developed for rock and soil mechanics to describe ductile fracture behavior over a wide range of loading states. Lou et al. [119] developed a micro-mechanism-inspired macroscopic ductile fracture criterion on the basis of assumptions that micro-voids proportionally nucleate with the equivalent plastic strain, the growth of voids was controlled by the stress triaxiality, and the coalescence of voids was governed by the maximum shear stress resulting in shear-linking up. In order to appropriately consider a reasonable cut-off value for the stress triaxiality, Lou et al. [120] proposed a cut-off value function which can be defined by Lode parameter, sensitivities to microstructures, temperature, strain rate, and so on. The uncoupled fracture criteria employ a scalar damage indicator, usually expressed in an integral form, which continually evolves with plastic deformation. In this case, the ductile fracture is known to initiate when the damage indicator reaches unity. By employing an uncoupled fracture criterion in the prediction of fracture initiation in aid of finite element analysis, Dunand and Mohr [63] made a great emphasis on the careful choice of hardening moduli in large plastic strain range in order to guarantee an accuracy of the fracture prediction. As an effort to evaluate the influence of material anisotropy on the fracture initiation, Luo et al. [123] introduced an uncoupled non-associated ductile fracture criterion based on the linear transformation of the strain tensor. Gu and Mohr [77] formulated an anisotropic extension of the Hosford–Coulomb shear localization criterion based on the linear transformation of the stress tensor. In modeling anisotropic ductile fracture criteria above, linear transformation tensors play a key role in the 4

construction of the anisotropic fracture configuration. The components of those transformation tensors were set as the model parameters to be calibrated. Various ductile fracture criteria developed are summarized in Tables 1 and 2.

1.2.2

Modeling of anisotropic yield behavior

For sheet metal forming and any other manufacturing processes involving plastic deformation, a mathematical description is essential to describe the stress state that will induce yielding or the onset of plastic deformation of the material. This phenomenological mathematical expression is called a yield criterion and many reviews concerning the development of yield criteria have been introduced in the literature [47, 143]. In developing constitutive laws for continuous media, following four basic principles are required [122]: (1) The principle of determinism: the behavior of materials at time t is determined by all the past history of the motion of all material points in the body until time t. (2) The principle of neighborhood : the behavior of a material point P at time t is determined by the behavior of an arbitrary small neighborhood. (3) The principle of coordinate invariance: the constitutive laws of materials are independent of the coordinate system (using tensor notation or abstract notation to meet naturally). (4) The principle of material objectivity: the constitutive laws of materials are independent of the observer who uses different clocks and does different movements (also known as frame non-difference principle). In other words, the constitutive laws are independent of the rigid motion of the spatial coordinates. The detailed discussion of the objectivity of stress and strain tensors can be found in Eringen [65,66], Guo [192], and Marsden and Hughes [126]. At an initial stage of the development of yield criteria, the simplest and most popular yield criteria were proposed to describe the isotropic material behavior: the Tresca [173]; the von Mises [177]; and the Hosford’s 72 non-quadratic yield criteria [91]. All these three yield criteria can be expressed in terms of the principal values of the Cauchy stress or the deviatoric stress as: f (sij ) = |S1 − S2 |m + |S2 − S3 |m + |S3 − S1 |m = 2¯ σm

(1.1)

For an exponent m = 2 or m = 4, Equation 1.1 reduces to the von Mises criterion, while for m = 1 or m = ∞, it leads to the Tresca criterion. Equation 1.1 itself represents the form of the Hosford’s 72 criterion and the exponent of m is in general set to 6 and 8 for BCC and FCC materials, respectively, for a good approximation of yield loci computed using the Bishop–Hill 5

crystal plasticity model [82, 117]. For sheet metals manufactured by cold rolling, the texture of the material tends to be oriented, which results in the directional material behavior. Generally, sheet metals are characterized by three mutually orthogonal planes of symmetry. These three directions are expressed as x, y, and z which correspond to the rolling, transverse and normal directions of the sheet, respectively. Taking into consideration of the orthotropic symmetry of the sheet metal, the following necessary and sufficient conditions should be satisfied in the formulation of an anisotropic yield criterion [157]: f (σxx , σyy , σzz , −σxy , σxz , −σyz ) = f (σxx , σyy , σzz , σxy , σxz , σyz )

(1.2a)

f (σxx , σyy , σzz , −σxy , −σxz , σyz ) = f (σxx , σyy , σzz , σxy , σxz , σyz )

(1.2b)

for any stress state. The first anisotropic yield criterion is introduced by Hill [83], which is called the Hill’s 48 criterion. The Hill’s 48 criterion is an extension of the von Mises criterion to orthotropic materials and is well suitable for special materials and textures especially for steel sheets; however, an anomalous behavior was reported by Pearce [144] and Woodthorpe and Pearce [181] for aluminum alloy sheets that gave r-values smaller than the unity. The anomalous behavior is a yield strength in the balanced biaxial tension tests that is greater than the uniaxial yield stress, which is contrary to the prediction of the Hill’s 48 criterion. Hill [84] proposed a non-quadratic yield criterion, the Hill’s 79 yield criterion, to describe materials other than steel and derived four special cases from the general form. All four cases encompass anomalous behavior for sheet metal having rn < 1 when proper values are assigned to the exponent of the Hill’s 79 criterion. The convexity of Hill’s 79 criterion was studied by Zhu et al. [194] who found that the predicted yield surfaces are concave or even unbounded for some combination of the exponent and r-values. The most widely used expression of the Hill’s 79 criterion is the case IV which applies to materials exhibiting planar isotropy for the plane stress state. Two other anisotropic yield criteria for the plane stress state were proposed by Hill [85] to consider the shear stress and the anomalous behavior. Independently from Hill, Hosford [90] generalized the Hosford’s 72 criterion for materials exhibiting orthotropic symmetry, which is called the Hosford’s 79 criterion. An important drawback of both the Hosford’s 79 and the Hill’s 79 criteria is that they do not consider the shear stress. Therefore, these yield criteria are not able to account for the continuous variation of plastic properties between the rolling and transverse directions of the sheet. Barlat and Lian [22] extended the Hosford’s 79 criterion so as to consider the influence of the shear stress, which is called the Yld89 or Barlat’s 89 criterion; however, the Yld89 criterion cannot describe the anomalous behavior. A six-component yield criterion was proposed by 6

Barlat et al. [21], which is called the Yld91 or Barlat’s 91 criterion, in order to accurately predict the yield stress in the balanced biaxial tension. There are only 4 anisotropic parameters in the Yld91 yield criterion which are generally evaluated by the r-values in the rolling, diagonal and transverse directions and the balanced biaxial tensile yield stress. In order to consider both yield stress and r-value directionalities, Barlat et al. [23] proposed another yield criterion, named as the Yld96 or Barlat’s 96 criterion, with 7 anisotropic parameters. The Yld96 criterion is able to describe the anisotropic of both the yield stress and the r-value. There are, however, several problems in the Yld96 criterion [20]: (1) There is no proof of convexity, which is an important requirement in numerical simulations to ensure the uniqueness of a solution. (2) The derivatives are difficult to obtain analytically, which is inconvenient for FE simulation. (3) The plane stress implementation in FE codes does not provide any particular problems and leads to good simulation results; however, some numerical problems for full stress states, which might be difficult to solve because of the relative complexity of the Yld96 formulation, have been encountered. To overcome these problems of the Yld96 criterion, Barlat et al. [20] newly developed another yield criterion, named as the Yld2000-2d criterion, which can guarantee convexity of the yield surface and make FE implementation and application simpler. There are 8 anisotropic parameters in the Yld2000-2d criterion, which are evaluated by σ0 , σ45 , σ90 , σb , r0 , r45 , r90 , and rb . Karafillis and Boyce [95] introduced the so-called Isotropic Plasticity Equivalent (IPE) theory with a more general yield criterion and a linear transformation that can accommodate other material symmetries. Based on the Yld89 criterion, the BBC2000 criterion was proposed by Banabic et al. [14], which has the same number of anisotropic parameters as that in the Yld2000-2D criterion. Most of the yield criteria reviewed above were, however, developed for the plane stress state. To describe the material anisotropy under the full stress state, Barlat et al. [18] developed two 3D yield criteria with 18 and 13 anisotropic parameters, respectively, based on the linear transformation of the Cauchy stress tensor. Recently, Cazacu et al. [42] developed a macroscopic orthotropic yield criterion that can describe both the anisotropy of a material and the yielding asymmetry between the tension and the compression.

7

1.3

Objective and scope of the thesis The objective of the thesis is to develop new anisotropic ductile fracture criterion capable

of describing the anisotropy effect on the equivalent plastic strain at the onset of fracture over a wide range of loading conditions as well as to suggest a single representative fracture forming limit diagram as a tool for evaluating the forming severity of an anisotropic material. Three kinds of fracture forming limit diagrams (FFLDs), such as a strain-based FFLD, a stressbased FFLD, and a polar effective plastic strain-based FFLD, will be reviewed and graphically represented based on a theoretical derivation from the anisotropic ductile fracture criterion to monitor the ductile fracture in AHSS sheets during sheet metal forming processes. A scaling method for construction of a single facture-based forming limit curve will be discussed to predict the onset of fracture for the anisotropic material. The anisotropy effect on the equivalent plastic strain at the onset of fracture will be evaluated for the DP980 1.2t (DP980 with the thickness of 1.2 mm) steel sheet from uniaxial tensile tests of various types of specimens including the pure shear, the uniaxial tension, and the plane strain tension. The fracture envelopes will be constructed according to the three different loading directions: the Rolling Direction (RD); the Diagonal Direction (DD); and the Transverse Direction (TD). A damage model will be introduced to calibrate the model parameters of the anisotropic ductile fracture criterion to consider the change of loading path in complicated forming processes. Damage evolutions of all loading cases will be computed using the calibrated parameters and strain histories utilized for construction of the anisotropic ductile fracture criterion to confirm whether the calibrated parameters are properly determined or not. Three-dimensional fracture envelopes of the DP980 1.2t steel sheet according to the three loading directions of 0°(RD), 45°(DD), and 90°(TD) will be represented in the space of (ηH , LP , ε¯f ) and (ηH , LP , σ ¯f ) based on the Hill’s 48 criterion, respectively. The corresponding fracture loci will be also extracted from the three-dimensional fracture envelopes to quantitatively predict the fracture strain in sheet metal forming processes. The FFLDs proposed will be reviewed with respect to the maximum principal stress direction to assess the forming severity and thus a representative single fracture-based forming limit criterion will be suggested to predict the onset of fracture for the anisotropic material. The developed fracture criterion will be applied to the numerical simulation by implementing it into the commercial Finite Element (FE) software of ABAQUS/Explicit by means of a Vectorized User-defined MATerial subroutine (VUMAT). Finally, the validity of the developed fracture criterion will be verified from comparison of the actual square cup drawing test and the FE analysis.

8

Chapter 2 Theoretical Backgrounds 2.1

Characterization of the material anisotropy For a more accurate evaluation of the sheet metal behavior during forming processes, it is

essential to acquire a quantitative description of the plastic material behavior by a constitutive model in numerical simulation. The von Mises criterion [177] is one of the simple and popular yield criteria immensely utilized for the analysis of plastic behavior for ductile materials, which is proposed based on the observation that hydrostatic pressure does not bring about plastic straining of the material. This criterion assumes that the transition from an elastic to a plastic state take places, independent of the stress state type when the elastic energy of distortion reaches a critical value. In 1948, Hill [83] proposed an anisotropic yield criterion as an extension of the von Mises criterion supposing that the material has an anisotropy with three orthogonal symmetry planes. Due to its simplicity, the Hill’s 48 criterion is widely adopted to describe the anisotropic deformation behavior of sheet metals. With increasing concerns about the accuracy of the prediction for anisotropic material behavior, dozens of complicated anisotropic yield criteria have been proposed such as those of the Hill’s family [83–86, 92, 108, 114], the Hershey’s family [13–15, 19, 21–23, 32, 52, 90, 95], the Drucker’s family [40–42], the polynomial family [51, 75, 157], and other yield criteria [67, 175]. The experimental quantities required for calibration of those yield criteria are briefly summarized in Table 3. Among them, the Hill’s 48 [83] criterion was chosen to describe the plastic anisotropy of DP980 1.2t steel sheet and to develop an anisotropic ductile fracture criterion for Advanced High-Strength Steel (AHSS) sheets due to the simplicity of its numerical formula and the small number of experimental data required for calibration. In this section, theoretical backgrounds of the Hill’s 48 criterion will be briefly reviewed.

2.1.1

Hill’s 48 criterion

Hill extended the von Mises criterion [177] to consider the plastic behavior of anisotropic polycrystalline metals on the basis of following assumptions [83]:

(1) The material is orthotropic; That is, there exist three orthogonal symmetry planes at each material point. The intersections of three planes are the principal axes of anisotropy;

9

(2) Hydrostatic pressure does not cause yielding; (3) There is no Bauschinger effect.

The Hill’s 48 criterion is expressed by a quadratic function of the following form according to the principal anisotropic axes of x, y, and z which usually coincide with the rolling, transverse and normal directions of the sheet: 2 2 2 2f (σij )≡F (σyy − σzz )2 + G(σzz − σxx )2 + H(σxx − σyy )2 + 2Lσyz + 2M σzx + 2N σxy = 1 (2.1)

where f is the yield criterion, F , G, H, L, M , and N are the parameters that characterize the current state of anisotropic yielding behavior. The material parameters of the Hill’s 48 criterion can be expressed as dependencies of the tensile yield stresses (X, Y, Z) and the shear yield stresses (R, S, T ) associated to the principal anisotropic directions: 1 F = 2



1 G= 2



1 H= 2



2L =

1 , R2

1 1 1 + 2− 2 2 Y Z X



1 1 1 + 2− 2 2 Z X Y



1 1 1 + 2− 2 2 X Y Z



2M =

1 , S2

(2.2a)

(2.2b)

2N =

(2.2c)

1 T2

(2.2d)

If L = M = N = 3F = 3G = 3H, Equation 2.1 reduces to the von Mises criterion for isotropic materials. In order to give a complete description of the anisotropic behavior of the material, six independent yield stresses, as well as the orientation of the principal anisotropic axes, need to be known. Under the plane stress condition (σzz = σxz = σyz = 0), Hill’s 48 criterion reduces to 2 2 2 2f (σij )≡(G + H)σxx − 2Hσxx σyy + (H + F )σyy + 2N σxy =1

(2.3)

In numerical simulation of the sheet metal forming process, the material parameters are denoted as r0 , r45 , r90 , and yield stresses in the directions of the principal anisotropic axes as X = σ0 , Y = σ90 . The relation among r-values and the material parameters F , G, H, and N can be easily obtained by the associative flow rule as follows: r0 =

H , G

r90 =

H , F 10

r45 =

N 1 − F +G 2

(2.4)

The material constants F , G, H, and N can be expressed in terms of σ0 , r0 , r45 , and r90 through F =

σ02 (1

G=

H=

N=

r0 + r0 )r90

(2.5a)

1 σ02 (1 + r0 )

(2.5b)

r0 + r0 )

(2.5c)

σ02 (1

(2r45 + 1)(r0 + r90 ) 2σ02 (1 + r0 )r90

(2.5d)

Equation 2.5 can be rewritten in terms of σ0 , σ45 , σ90 , and r0 as follows: F =

1 r0 − 2 2 σ90 σ0 (1 + r0 )

G=

H=

N=

(2.6a)

1 + r0 )

(2.6b)

r0 + r0 )

(2.6c)

σ02 (1

σ02 (1

2 1 (1 − r0 ) 2 − 2σ 2 (1 + r ) − 2σ 2 σ45 0 0 90

(2.6d)

According to Equation 2.5, Equation 2.3 can be re-written as 2 σxx −

r0 (1 + r90 ) 2 r0 + r90 2r0 2 σxx σyy + σ + (2r45 + 1)σxy = σ02 1 + r0 r90 (1 + r0 ) yy r90 (1 + r0 )

(2.7)

In the case that principal axes of the yield criterion coincide with the principal anisotropic axes, i.e. σxx = σ1 , σyy = σ2 , Equation 2.7 can be expressed as σ12 −

2r0 r0 + r90 2 σ1 σ2 + σ = σ02 1 + r0 r90 (1 + r0 ) 2

(2.8)

From Equation 2.8, it can be concluded that four experimental results, e.g. σ0 , r0 , r45 , and r90 , are needed to calibrate the yield criterion under the plane stress condition. For a normal anisotropic material (r0 = r45 = r90 = rn ), Equation 2.8 reduces to 11

σ12 −

2rn σ1 σ2 + σ22 = σu2 1 + rn

(2.9)

where σu is the uniaxial yield stress. In the case that if rn < 1, the yield locus predicted by the Hill’s 48 criterion is located inside the one given by von Mises yield criterion; if rn > 1, the Hill’s 48 yield locus is outside the von Mises yield locus. For the equi-biaxial tension case (σ1 = σ2 = σb ), the relation between the uniaxial and equi-biaxial yield stresses can be obtained by rewriting Equation 2.9 through r σb = σu

1 + rn 2

(2.10)

where σb is the equi-biaxial yield stress. It is obvious that if rn > 1, then σb > σ0 and if rn < 1, then σb < σ0 . However, Pearce [144] and Woodthrope and Pearce [181] observed that some materials, in particular aluminum alloys, show σb > σ0 even though rn < 1. This has been called ‘anomalous behaviour’ and the anisotropic behaviour for materials exhibiting this anomalous behaviour cannot be properly described by Hill’s 48 criterion. Moreover, the Hill’s 48 criterion cannot represent the particular case of r0 /r90 > 1 when σ0 /σ90 < 1 (or vice-versa) because the criterion predicts σ0 = σ90

s

r0 (1 + r90 ) r90 (1 + r0 )

(2.11)

This has been called the ‘second order anomalous behavior’ [16]. The Hill’s 48 criterion can only be applied to materials forming four ears in axisymmetric deep drawing processes whereas in practice also different numbers of ears are observed [11]. Except for some particular cases, the Hill’s 48 criterion is widely utilized in numerical simulations of auto-body steel sheets due to the simplicity of its numerical formula and it has been known that Hill’s 48 criterion shows a good performance to describe the anisotropy of most steel sheets.

2.2

Anisotropic coefficients calculation To determine the anisotropic characteristic parameters of yield criteria, several experimen-

tal data are required such as the yield stresses and the r-values evaluated from different loading directions. In the case of the Hill’s 48 criterion, all the anisotropic parameters can be computed explicitly using Equation 2.5 whereas some advanced yield criteria need numerical iteration to determine the anisotropic parameters. For example, the anisotropic parameters of the Yld20002d criterion cannot be computed explicitly and have to be obtained through numerical methods due to the complexity of its numerical formula. The procedure to calculate the anisotropic parameters of the Hill’s 48 criterion will be explained in detail in the following section. 12

2.2.1

Calculation of the yield stress in an arbitrary direction

Considering that a uniaxial tensile test is carried out with a specimen cut at the angle θ to the rolling direction, the stress components in the orthotropic system can be expressed as σxx = σθ cos2 θ

(2.12a)

σyy = σθ sin θ2

(2.12b)

σxy = σθ sin θ cos θ

(2.12c)

where σθ denotes the yield stress in the tensile direction. The yield criterion given in Equation 2.3 then becomes   2f (σij ) = (G + H) cos4 θ − 2Hsin2 θcos2 θ + (F + H) sin4 θ + 2N sin2 θcos2 θ σθ2 = 1

(2.13)

Solving for σθ in the above equation leads to  σθ =

1 2 2 4 (G + H) cos θ − 2Hsin θcos θ + (F + H) sin4 θ + 2N sin2 θcos2 θ

1 2

σ ¯H

(2.14)

Applying θ = 0°, θ = 45°, and θ = 90° into Equation 2.14, the relation between σ0 , σ45 , σ90 , and the anisotropic parameters of the Hill’s 48 criterion is obtained as follows: σ ¯H

σ0 =

(2.15a)

1

(G + H) 2 2¯ σH

σ45 =

1

(2.15b)

(G + F + 2N ) 2 σ ¯H

σ90 =

1

(2.15c)

(F + H) 2

2.2.2

Calculation of the r-value in an arbitrary direction

The r-value, called Lankford parameter, of a specimen subjected to tensile loading in θdirection from the rolling direction is defined by rθ =

ε˙θ+90° ε˙θ+90° ε˙θ+90° =− =− ε˙t ε˙θ + ε˙θ+90° ε˙xx + ε˙yy

(2.16)

Using the associative flow rule, ∂f ε˙θ = λ˙ ∂σθ 13

(2.17)

and noting that ∂σxx = cos2 (θ + 90°) = sin2 θ ∂σθ+90°

(2.18a)

∂σyy = sin2 (θ + 90°) = cos2 θ ∂σθ+90°

(2.18b)

∂σxy = sin(θ + 90°) cos(θ + 90°) = − sin θ cos θ ∂σθ+90°

(2.18c)

We, then, can, by means of the chain rule of differentiation, express rθ as ∂f ∂f ∂f sin2 θ + cos2 θ − sin θ cos θ ∂σxx ∂σyy ∂σxy rθ = − ∂f ∂f + ∂σxx ∂σyy

(2.19)

where differentiations of the yield criterion are ∂f = G (σxx − σyy ) + Gσxx ∂σxx

(2.20a)

∂f = F σyy + H (σyy − σxx ) ∂σyy

(2.20b)

∂f = −Gσxx − F σyy ∂σzz

(2.20c)

∂f = 2N σxy ∂σxy

(2.20d)

Substituting Equation 2.20 into Equation 2.19 gives rθ =

H + Htan2 θsin2 θ + 2N sin2 θ − (F + G + 3H) sin2 θ G + F tan2 θ

(2.21)

Applying θ = 0°, θ = 45°, and θ = 90° into Equation 2.21, the relation between the r-values of r0 , r45 , and r90 and the anisotropic parameters of the Hill’s 48 criterion is obtained as follows: H G

(2.22a)

2N − (F + G) 2(F + G)

(2.22b)

r0 =

r45 =

14

r90 =

H F

(2.22c)

It is worth to mention that the Hill’s 48 criterion may not be always convex, and thus it is necessary to confirm whether the values of the anisotropic parameters satisfy the inequalities given in Equation 2.23 to ensure that the yield criterion is real and convex [83]: F + G > 0,

N > 0,

F G + GH + HF > 0

15

(2.23)

Chapter 3 Quantification of Plastic Anisotropy at Various Orientations 3.1

Tensile tests at various orientations In order to investigate the anisotropy effect on the equivalent plastic strain at the onset

of fracture especially for AHSS sheets, the DP980 1.2t steel sheet was tested and analyzed in the present paper. The chemical compositions of the DP980 1.2t steel sheet tested are given in Table 4. DP(Dual-Phase) steels usually consist of a ferrite matrix(soft phase, offers ductility) containing a martensite second phase (hard phase, offers strength) in the form of islands. DP steels can be easily deformed due to its low yield to tensile strength ratios of around 0.50–0.60, and it exhibits a high level of total elongation close to that of TRIP(TRansformation Induced Plasticity) steel. From these characteristics of the DP steel, cold rolled DP sheets are generally utilized for outer panels, structural parts, and crashworthiness parts of an auto-body. Three different types of tensile tests were carried out to quantitatively measure the equivalent plastic strain at the onset of fracture of the DP980 1.2t steel sheet over a wide range of stress states. Tensile specimens were designed to have various shapes, as shown in Figure 4, to induce the typical loading conditions at the material point where the fracture initiates from uniaxial tensile tests. The three types of tensile specimens are: classical dog-bone; pure shear; and flat-grooved plane strain specimens. Each specimen was fabricated along the three different orientations: the rolling direction; the diagonal direction; the transverse direction of sheet metal to evaluate the anisotropy effect on the equivalent plastic strain at the onset of fracture. The universal testing machine of INSTRON 5583, as shown in Figure 5, was utilized for the quasistatic tensile tests. Deformation images on the surface of specimens were acquired by means of a high-speed camera of Phantom v9.0 to measure the strain distribution and strain path history. Various mechanical properties of the DP980 1.2t steel sheet were obtained from uniaxial tensile tests at the three loading directions: Yield Stress (YS); Ultimate Tensile Strength (UTS); uniform elongation (Elu ); and fracture elongation (Elf ). Yield stresses were measured with the 0.2 % offset method from the obtained engineering stress–engineering strain curves. The hardening behaviors at the three loading directions are represented in Figure 6 with the curves interpolated by a mixed Swift–Voce strain hardening model as given in Equation 3.1.

16

 p  σ ¯ = α K(ε0 + ε¯p )n + (1 − α) A + B 1 − e−C ε¯ | {z } {z } | Swift law [169]

(3.1)

Voce law [176]

The mechanical properties of the DP980 1.2t steel sheet obtained from uniaxial tensile tests and the parameters of the mixed Swift–Voce hardening model obtained are summarized in Tables 5 and 6, respectively. Equivalent plastic strains at the onset of fracture and strain path histories were measured at the material point which has a maximum equivalent plastic strain just before the fracture initiation to calibrate the parameters of fracture criterion with the aid of the two-dimensional digital image correlation (DIC) method. Here, the DIC method is now extensively being utilized to accurately measure the surface deformations and stresses [62, 64, 71, 72, 102, 123, 152, 171] since it was first applied to experimental mechanics by [48], Ranson and Peters [148], and Peters et al. [145]. As a preparation for applying the DIC method, a random speckle pattern with white and black spray-paints was generated on the surface of the specimen to enhance the contrast of captured images. The optical measurement system was set up to acquire the deformation images of the specimen with the FASTCAM SA4 motion analysis camera during the uniaxial tensile tests. The grayscale digital images captured were analyzed by the DIC method in aid of a commercial digital image processing software, ARAMIS v6.3.0. The spatial resolution and frame rate used were around 0.02 mm/pixel and 20 frames/sec, respectively. With adopting this methodology for strain measurement, tensile tests with three types of specimens were carried out three times for each type to confirm the reproducibility. The equivalent plastic strains at the onset of fracture were evaluated from tensile tests as shown in Figures 7 and 8, respectively. In comparison of tests results, representative strain paths for each loading direction were chosen from three test results, which satisfies a consistent strain path desired during tensile tests. The equivalent plastic strain at the onset of fracture was determined with the one close to the average value of three test results for typical loading condition. Figure 9 shows the representative strain paths determined and Table 7 explains the quantitative values of the equivalent plastic strain at the onset of fracture evaluated from each test.

3.2

Mesurement of the r-value at various orientations For sheet metals, the variation of plastic behavior with direction and deep-drawability are

usually assessed by means of r-value which is also called Lankford parameter or anisotropy coefficient [107]. The r-value is a quantity describing the ability of a sheet metal to resist thinning or thickening when subjected to either tensile or compressive load in the plane of the

17

sheet as explained in ASTM (American Society for Testing and Materials) E517 standard [2]. This parameter is determined by uniaxial tensile tests on sheet specimens. The r-value is defined as r=

ε˙22 ε˙33

(3.2)

where ε˙22 and ε˙33 are the strain increments along the width and thickness directions, respectively. This parameter should not be confused with the Poisson’s ratio which is defined as a ratio of the transverse contraction strain to the longitudinal extension strain in the direction of stretching force (ν = −ε˙22 /ε˙11 ). In the case that a material is fully isotropic, the r-value becomes a unity, which indicates that the width and thickness strains have the same value. If the r-value is greater than unity, the width strains will be dominant resulting in an increase in thinning resistance, such that drawing of the material from the blank is easier than reducing the thickness of the sheet. On the other hand, for the material exhibiting the r-value less than unity, the thickness reduction will dominate in stamping and drawing operations. In the case that a material has planar anisotropy, the r-value depends on the in-plane direction of the sheet metal. Considering a specimen cut at the angle θ to the rolling direction, the instantaneous r-value rθ is defined as a ratio of the plastic strain rates associated with the width direction ε˙θ+90 and the thickness direction ε˙33 [11]: rθ =

ε˙θ+90° ε˙33

(3.3)

Due to practical difficulties in measuring the through-thickness plastic strain in sheet metals, the r-value is conventionally determined with the postulate of the plastic incompressibility of metals, i.e. ε˙ii = 0. The r-value can be, hence, obtained from the measurement of in-plane strain rates as follows: rθ = −

(ε˙θ+90° /ε˙θ ) ε˙θ+90° =− ε˙θ+90° + ε˙θ 1 + ε˙θ+90° /ε˙θ

(3.4)

where ε˙θ denotes the longitudinal strain increment at the angle θ to the rolling direction. The normal and planar anisotropy coefficients can be then defined as rn =

r0 + 2r45 + r90 r0 − 2r45 + r90 , ∆R = 4 2

18

(3.5)

The normal and planar anisotropy coefficients are directly related to the deep-drawability and the earing generation of the material. A high normal anisotropy coefficient increases the in-plane formability of the material whereas a high planar anisotropy coefficient increases the height of ears in a drawn steel cup. Typical anisotropy coefficients for various materials are summarized in Table 8 [76]. In the present paper, r-values of the DP980 1.2t steel sheet were measured in the range of uniform elongation by means of an optical system incorporated with the 2D-DIC method using specimens cut at the angle of 0°, 45°, and 90° to the rolling direction of the sheet metal. Test results are represented in Figure 10 and the r-values calculated are given in Table 5 in an average sense. The instantaneous r-value is also evaluated to investigate the evolution of anisotropy with the plastic deformation. To predict the evolution of the r-values appropriately, the instantaneous r-values of each loading direction are fitted by an asymptotic regression model which has a form of rθ (εθ ) = Aθ − Bθ exp (−Cθ εθ )

(3.6)

Figure 11 shows the instantaneous r-values fitted and their corresponding residuals. The parameters of the asymptotic regression model are summarized in Table 9 according to each loading direction. The transverse plastic strain–longitudinal plastic strain curves are re-obtained from the curves fitted and comparison of the obtained curves with the experimental data is conducted by calculating the residuals as shown in Figure 12. The relation between the transverse plastic strain and the longitudinal plastic strain can be obtained from Equation 3.4 as: ε˙θ+90° rθ =− → εθ+90° = ε˙θ 1 + rθ

Z

tu

Z ε˙θ+90° dt = −

0

0

tu

rθ ε˙θ dt 1 + rθ

(3.7)

where tu denotes the time elapsed during uniform elongation of the specimen. It is noted from the comparison that the fitted curves predict the material behaviors with great accuracy, showing that the absolute values of the residuals for each loading direction are less than about 0.0002. There shows little difference between the r-values obtained directly by the transverse plastic strain–longitudinal plastic strain curves in an average sense and the ones obtained after saturation of the instantaneous r-values which are fitted by the asymptotic regression model. The saturated r-values are listed in Table 5. Figure 13 represents the instantaneous normal and planar anisotropic quantities evaluated using the instantaneous r-values having the same plastic work. With the use of r-values evaluated above, the anisotropic characteristic parameters of the Hill’s 48 criterion were calibrated and the parameters obtained are summarized in Tables 10 and 11. In calibrating the Hill’s 48 anisotropic parameters, the material anisotropy in the off-axes 19

directions was neglected. Figure 14 represents the Hill’s 48 yield loci with the von Mises yield locus to confirm the level of anisotropy.

20

Chapter 4 Characterization of the Stress State 4.1

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηv , θ, σ ¯v ) Cauchy stress tensor σ has six components at a material point on the basis of symmetry

conditions from the balance of angular momentum, and thus the material state can be interpreted in six-dimensional space when subjected to certain loading conditions. It is, however, quite complicated to directly give a comprehensive understanding of the stress states that the material undergoes. This can be released by dealing with the material state using three principal stress components of (σ1 , σ2 , σ3 ) in Haigh–Westergaard space where the coordinates are three principal stresses, such that a specific stress state is represented with a specific stress vector as shown in Figure 15(a). In this space, the orientation of the principal stress components is not considered. In Haigh–Westergaard space, the plane, which is passing through the origin O and is −→ perpendicular to a hydrostatic axis, is called a π-plane. The stress vector OP can be decomposed −−→ −−→ into two components as hydrostatic part OO0 and deviatoric part O0 P. The vector component −−→0 −−→ OO is perpendicular to the π-plane and the vector O0 P is on the π-plane, respectively. The −−→ magnitude of the vector OO0 is expressed in terms of the hydrostatic stress as √ −−→ |OO0 | = 3σm

(4.1)

where σm is the mean stress and is defined through σm =

1 (σ1 + σ2 + σ3 ) 3

−−→ The magnitude of the vector O0 P is linearly related to the equivalent stress; i.e., r −→ 2 σ ¯v |OP| = 3

(4.2)

(4.3)

where σ ¯v is the von Mises equivalent stress and is expressed as 1 σ ¯v = √ 2

q (σ1 − σ2 )2 + (σ2 − σ3 )2 + (σ3 − σ1 )2

(4.4)

−→ When the principal stress vector OP is projected on the deviatoric plane as shown in Figure 15(b), an angle θ can be defined as an indicator, called a Lode angle, differentiating the stress 21

state between tension, shear, and compression. The influence of the Lode angle is usually used as Lode dependence in literature and there exist different definitions for this indicator. The Lode angle is defined as the angle θa , as shown in Figure 15(b), to the positive direction of the σ1 = −σ3 axis [57, 70, 183]. On the one hand, the Lode angle is expressed as the positive angle to the positive direction of σ1 axis [8, 46, 137]. For the transformation of the principal stress state of (σ1 , σ2 , σ3 ) into the space of (σm , θ, σ ¯v ) which is related to the stress invariants, the latter one is adopted in this study. The Lode angle is defined by a relative ratio of principal deviatoric stresses. Each definition of the Lode angle can be expressed as −1

θa = tan



 1 √ (2χ − 1) , 3

−1

θ = cot



2 √ 3



1 1 − χ 2

 (4.5)

where the relative ratio of principal deviatoric stresses is defined as [182, 188, 189] χ=

σ2 − σ3 σ 0 − σ30 = 20 σ1 − σ3 σ1 − σ30

(4.6)

With the quantities (σm , θ, σ ¯v ) discussed above, the cylindrical dimension can be described in Haigh–Westergaard space to indicate the principal stress state as follows: 2 ¯v cos θ σ1 = σm + σ10 = σm + σ 3   2 2 0 ¯v cos π−θ σ2 = σm + σ2 = σm + σ 3 3   2 4 0 σ3 = σm + σ3 = σm + σ ¯v cos π−θ 3 3

(4.7a) (4.7b) (4.7c)

The three cylindrical dimensions of σm , σ ¯v , and θ can be expressed in terms of stress invariants through 1 1 σm = I1 = σkk 3 3

(4.8)

r p σ ¯v = 3J2 = 1 θ = cos−1 3

3 0 0 σ σ 2 ij ij ! √ 3 3J3 3/2 2J2

where

0 σij = σij − σm δij

where

1 0 0 0 J3 = σij σjk σkl 3

(4.9)

(4.10)

where I1 is the first invariant of the Cauchy stress tensor σ, J2 and J3 are the second and the invariants of the stress deviator tensor σ 0 , respectively. With these stress invariants, the stress state at a material point can be uniquely defined.

22

In order to obtain the point of view of stress states with a direction of the principal stress vector, it is convenient to work with the hydrostatic stress normalized by the von Mises equivalent stress, called stress triaxiality, as in Equation 4.11 since it is linearly related to a geometric −→ angle between the hydrostatic axis and the principal stress vector OP in the cylindrical coordinate system: ηv =

σm σ ¯v

(4.11)

Through a combination of this dimensionless pressure ηv and the Lode angle, the direction of principal stress vector can be uniquely determined. As a result, it can simply say that the principal stress states lying in the same direction in Haigh–Westergaard space are regarded as having the same loading condition. The principal stress state now can be uniquely defined with use of the von Mises equivalent stress which defines the magnitude of the principal stress vector in Haigh–Westergaard space as follows: 

 2 σ1 = σm + s1 = ηv + cos θ σ ¯v 3    2 2 π−θ σ ¯v σ2 = σm + s2 = ηv + cos 3 3    2 4 σ3 = σm + s3 = ηv + cos π−θ σ ¯v 3 3

(4.12a) (4.12b) (4.12c)

The Lode angle θ has a relation with the normalized third deviatoric invariant ξ through √ 3 3 J3 27 J3 = ξ= = cos 3θ 2 σ ¯v3 2 J 3/2 2

(4.13)

For a more detailed derivation, readers are recommended to refer to Khan et al. [97]. The normalized third deviatoric invariant, hereinafter, will be denoted as a Lode angle parameter, as it is a function of the Lode angle θ. Since the range of the Lode angle θ is 0 ≤ θ ≤ π/3, that of the Lode angle parameter ξ and the relative ratio of the principal deviatoric stresses χ are −1 ≤ ξ ≤ 1 and 0 ≤ χ ≤ 1, respectively. It can be simply showed that ξ = 1 (θ = 0, χ = 0) corresponds to axisymmetric tension or equi-biaxial compression, ξ = 0 (θ = π/6, χ = 0.5) corresponds to generalized shear condition (plane strain). The lower limit value ξ = −1, (θ = π/3, χ = 1) corresponds to axisymmetric compression or equi-biaxial tension as shown in Figure 15(b). In some research [9, 124], the Lode angle parameter is alternatively defined as 2 6θ = 1 − cos−1 ξ θ¯ = 1 − π π

23

(4.14)

¯ have the same values for It is noted that the two Lode angle parameter definitions (ξ, θ) the limiting bounds of each typical stress states discussed above whereas the definitions for the midrange slightly differ from each other. The stress triaxiality has a relation with the stress invariants of I1 and J2 , whereas the Lode angle parameter ξ is formulated with the stress invariants of J2 and J3 . Under the proportional loading condition, the stress triaxiality and the Lode angle parameter remain constant during loading since those are related to invariants denoting the typical stress state. Wierzbicki and Xue [180] showed that there is a unique relation between the stress triaxiality and Lode angle parameter under the plane stress condition (σ3 = 0) through   27 1 2 ξ = − ηv ηv − 2 3

(4.15)

In industrial applications especially for sheet metal forming, the formulation under the plane stress condition is generally used, which requires only one of two indicators (ηv , ξ) to determine the stress state at a material point. The influence of these indicators on the material behavior gives a significant physical viewpoint in developing a fracture criterion because the stress triaxiality and the Lode angle parameter describe the pressure dependence and the Lode dependence on the material yielding, respectively.

4.2

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηv , LP , σ ¯v ) As discussed in section 4.1, the Cauchy stress state can be explained in Haigh–Westergaard

space (σ1 , σ2 , σ3 ). The typical principal stress state can be interpreted with a new indicator, called Lode parameter [116], which is introduced to represent the status of the intermediate principal stress as follows: LP =

2σ2 − σ1 − σ3 =2 σ1 − σ3



σ2 − σ3 σ1 − σ3



 − 1 ← σ2 =

σ1 + σ3 2



 + LP

σ1 − σ3 2

 (4.16)

where the principal stresses are in the order of σ1 ≥ σ2 ≥ σ3 . From relations between the three indicators (ηv , LP , σ ¯v ), the principal stress state can be uniquely defined since the Lode parameter is linearly related to the Lode angle θ. Using Equations 4.4 and 4.16, the following relation is obtained: σ1 − σ3 σ ¯v =q 2 L2P + 3 24

(4.17)

The simple mathematical manipulation through Equations 4.4, 4.16, and 4.17 gives   σ1 + σ3 LP L P σ = σm − q σ ¯v = ηv − q ¯v 2 2 2 3 LP + 3 3 LP + 3

(4.18)

From Equations 4.17 and 4.18 with a condition of the first deviatoric stress invariant J1 = 0, the principal stress state is uniquely defined in the space of (ηv , LP , σ ¯v ) as follows: 



3 − LP  3 − LP σ ¯v = ηv + q σ ¯v σ1 = σm + s1 = σm + q 3 L2P + 3 3 L2P + 3

(4.19a)

  2LP 2L P σ σ2 = σm + s2 = σm + q σ ¯v = ηv + q ¯v 2 2 3 LP + 3 3 LP + 3

(4.19b)

  3 + L 3 + LP P  σ ¯v σ ¯v = ηv − q σ3 = σm + s2 = σm − q 2 2 3 LP + 3 3 LP + 3

(4.19c)

For the plane stress condition, the relation between the stress triaxiality and the Lode parameter is uniquely defined from Equation 4.19. With the associative flow rule, the strain paths for isotropic materials can be solely expressed in terms of the Lode parameter through ε˙p1 : ε˙p2 : ε˙p3 = 2σ1 − σ2 − σ3 : 2σ2 − σ3 − σ1 : 2σ3 − σ1 − σ2 = 3 − LP : 2LP : −3 − LP

(4.20)

Subsequently, the principal strain rates, which are normalized by the von Mises equivalent strain rate, can be defined as ε˙p1 3 − LP = q p ε¯˙ 2 L2P + 3

(4.21a)

2LP ε˙p2 = q ˙ε¯p 2 L2P + 3

(4.21b)

ε˙p3 −3 − LP = q p ε¯˙ 2 L2P + 3

(4.21c)

The normalized principal strain rates in Equation 4.21 strictly comply with the order of ε˙p1 ≥ ε˙p2 ≥ ε˙p3 . The Lode parameter then can also be defined by the plastic strain increments through LP =

2ε˙p2 − ε˙p1 − ε˙p3 3ε˙p2 = ε˙p1 − ε˙p3 ε˙p1 − ε˙p3 25

(4.22)

As discussed in sections 4.1 and 4.2, the principal stress state (σ1 , σ2 , σ3 ) can be characterized in the space of (ηv , LP , σ ¯v ) as well as (ηv , θ, σ ¯v ). From these transformations, we can simply say that there is a certain relation between the Lode parameter and the Lode angle. By putting Equation 4.19 into Equation 4.13, the relation can be derived as ξ = cos 3θ =

LP (LP − 3) (LP + 3) 3/2 L2P + 3

(4.23)

The relation between the Lode parameter and the normalized Lode angle θ¯ is subsequentely obtained by the combination of Equations 4.14 and 4.23 as sin

π θ¯ LP (LP − 3) (LP + 3) = 3/2 2 L2P + 3

(4.24)

Equation 4.23 can be simply expressed from a relation between the maximum principal stresses in Equations 4.12 and 4.19 through 3 − LP cos θ = q 2 L2P + 3

(4.25)

Reversely, the Lode parameter can be calculated by the Lode angle as below: √ 3 tan θ − 3 √ LP = tan θ + 3

(4.26)

From Equations 4.14 and 4.25, the relation between the normalized Lode angle and the Lode parameter can be described in a simple manner as sin

π θ¯ LP = −q 6 L2 + 3

(4.27)

P

Another Lode angle related parameter χ as in Equation 4.6 can be expressed in terms of the Lode parameter as below: χ=

1 + LP 2

(4.28)

The values of each indicator, corresponding to typical stress states, are summarized in Table 12. Figure 16(a) represents the dependence of the stress triaxiality on various deviatoric stress state variables under the plane stress condition. The relations between each indicator are graphically represented in Figure 16(b) according to the normalized third invariant of the 26

deviatoric stress tensor. It is worth to note that the mechanism of void growth is limited under the low stress triaxiality (ηv ≤ 0) where the shear localization plays a main role in the fracture initiation. Under this failure mechanism, the size of void can be much smaller than that for failure governed by the growth of void and inter ligament necking (ductile fracture).

4.3

Transformation of the of the stress state (σ1 , σ2 , σ3 ) into the space of (ηH , LP , σ ¯H ) ¯ the stress triaxiality ηv plays With the Lode parameter LP or the normalized Lode angle θ,

a key role in modeling uncoupled ductile fracture criterion for the reason that the combination of these stress invariants represents the direction of three-dimensional principal stress vector (σ1 , σ2 , σ3 ) in Haigh–Westergaard space. It is, however, challenging to consider the directionality of typical stress states on the material orientation using these stress invariants because they are isotropic indicators. In this section, a new indicator is proposed in an attempt to involve the influence of material anisotropy on the stress triaxiality. The new indicator has a similar form of the stress triaxiality (ηv = σm /¯ σv ) but is defined with use of the Hill’s 48 equivalent stress for the normalization of the mean stress. Hereinafter, the new indicator will be also called the stress triaxiality for the sake of brevity. The general Hill’s 48 equivalent stress σ ¯H and the stress triaxiality defined by the Hill’s 48 equivalent stress ηH are defined as follows: 2 2 2 2 σ ¯H = F (σyy − σzz )2 + G(σzz − σxx )2 + H(σxx − σyy )2 + 2Lσyz + 2M σzx + 2N σxy

ηH =

σxx + σyy + σzz σm 1 = p σ ¯H 3 F (σyy − σzz )2 + G(σzz − σxx )2 + H(σxx − σyy )2

(4.29) (4.30)

In aid of the tensor transformation using the Rodrigues’ rotation formula [151] given by Equation 4.31a, the Cauchy stress state, defined in the reference coordinate system along the three orientation of sheet metal (RD, TD, and ND), can be expressed in terms of principal stresses, i.e., 

cos θr + (1 − cos θr )x2

(1 − cos θr )xy + z sin θr (1 − cos θr )xz − y sin θr



    R = (1 − cos θr )xy − z sin θr cos θr + (1 − cos θr )y 2 (1 − cos θr )yz + x sin θr    (1 − cos θr )xz + y sin θr (1 − cos θr )yz − x sin θr cos θr + (1 − cos θr )z 2 



  [σij ] = [RiI ][σIJ ][RjJ ] = Dσ  Eσ 27

Dσ Eσ Bσ Fσ

(4.31a)



  Fσ   Cσ

(4.31b)

with   σxx σxy σxz     [σij ] = σyx σyy σyz  ,   σzx σzy σzz

  σ1 0 0     [σIJ ] =  0 σ2 0    0 0 σ3

(4.31c)

    2 + c)2 2 2   A (tx (txy + sz) (txz − sy)   σ             2 2 + c)2 2   B (txy − sz) (ty (tyz + sx)     σ           σ1    Cσ       (txz + sy)2 (tyz − sx)2 (tz 2 + c)2  σ =    2  (4.31d)   2 + c)(txy − sz) 2 + c)      D (tx (txy + sz)(ty (txz − sy)(tyz + sx) σ     σ      3       2 + c)(txz + sy) 2 + c)     E (tx (txy + sz)(tyz − sx) (txz − sy)(tz   σ          2 2 Fσ (txy − sz)(txz + sy) (ty + c)(tyz − sx) (tyz + sx)(tz + c) ˆ is a unit vector for a rotation axis where t = 1 − cos θr , c = cos θr , s = sin θr , u ˆ = xˆi + yˆ j + zk and θr is the rotation angle as shown in Figure 17. By substituting Equations 4.16 and 4.31b into Equation 4.30, the stress triaxiality defined by the Hill’s 48 equivalent stress can be stated in terms of the Lode parameter and principal stresses of σ1 and σ3 . With this transformation, the anisotropy effect of a certain loading state on the material orientation can be considered. The Hill’s 48 equivalent stress and the stress triaxiality are subsequently defined as Equation 4.32a and Equation 4.32b, respectively. 2 2 2 2 σ ¯H = F (σyy − σzz )2 + G(σzz − σxx )2 + H(σxx − σyy )2 + 2Lσyz + 2M σzx + 2N σxy 2

2

2

2

2

2

(4.32a)

2

= (σ1 − σ3 ) [F f + Gg + Hh + 2(Ll + M m + N n )] ηH =

σxx + σyy + σzz 1 σm = p σ ¯H 3 (σ1 − σ3 )2 [F f 2 + Gg 2 + Hh2 + 2(Ll2 + M m2 + N n2 )]

(4.32b)

with 

  1 − 2x2 + x4 − 2y 2 + c2 x4 + 2x2 y 2 + 2cx2 − 2cx4 − 4cx2 y 2 + 2c2 x2 y 2 − 2xyzs + 2xyzsc        1 − LP 2 2 2   + −c + 2c y − 2xyzs + 2xyzsc f =  (4.32c) 2        1 + LP + −x2 − y 2 + 2y 4 + c2 x2 − c2 y 2 + 2c2 y 4 + x2 y 2 + 4cy 2 − 4cy 4 + c2 x2 y 2 − 2cx2 y 2 2 

x2 − 2x4 + y 2 + c2 x2 − 2c2 x4 − c2 y 2 − x2 y 2 − 4cx2 + 4cx4 − c2 x2 y 2 + 2cx2 y 2





      1 + LP  −1 + 2x2 + 2y 2 − y 4 − c2 y 4 − 2x2 y 2 − 2cy 2 + 2cy 4 − 2c2 x2 y 2 + 4cx2 y 2 − 2xyzs + 2xyzsc  g = +  2        1 − LP + −c2 + 2c2 x2 + 2xyzs − 2xyzsc 2

  x4 − 2c2 x2 + c2 x4 + 2cx2 − 2cx4 + 2xyzs − 2xyzsc        1 + LP 4 2 2 2 4 2 4   + −y + 2c y − c y − 2cy + 2cy + 2xyzs − 2xyzsc h=  2      1 − LP 2 2 2 2 2 2 2 2 2 2 2 2 2 2 + −1 + 2c + x + y − c x − c y − x y − c x y + 2cx y 2

(4.32d)



28

(4.32e)

    l=   

−xs − yz + x3 s + x2 yz − cx3 s + xy 2 s − cxy 2 s − 2cx2 yz + c2 x2 yz    1 + LP + y 3 z − 2cy 3 z + c2 y 3 z + cyz 2    1 − LP + c2 yz + csx + xy 2 s − cxy 2 s 2



x3 z − 2cx3 z + c2 x3 z + cxz

       

(4.32f)





      3 s + xy 2 z + cy 3 s − x2 ys + cx2 ys − 2cxy 2 z + c2 xy 2 z  + 1 + LP ys − xz − y m=  (4.32g) 2        1 − LP 2 2 2 −c xz + cys + x ys − cx ys + 2 

  x3 y + cxy − c2 xy − 2cx3 y − x2 zs + c2 x3 y + cx2 zs     n =   1 + LP   (4.32h)  1 − LP 3 2 3 2 2 3 2 + xy + cxy − c xy − 2cxy + y zs + c xy − cy zs + (−csz) 2 2

For the application of the proposed stress triaxiality to sheet metal forming, the unit vector ˆ = k. ˆ This makes it possible for the stress of a rotation axis can be defined as u ˆ = xˆi + yˆ j + zk triaxiality to have an explicit viewpoint about the dependence of the typical stress states on the material orientation with the consideration of the maximum principal stress direction. Following the above procedure, Equation 4.32b can reduce to ηH =

σm 1 σ1 + σ2 + σ3 √ = σ ¯H 3 (σ1 − σ3 ) T

(4.33a)

where 

2



1 + LP + G cos θS + sin2 θS T =F 2    N 1 − LP 2 2 2 + H cos 2θS + sin 2θS 2 2 1 + LP cos2 θS + sin2 θS 2

2

2 (4.33b)

In this case, the rotation angle of θr is denoted as θs which represents the the maximum principal stress direction to the reference rolling direction as shown in Figure 75. After the simple mathematical manipulation with Equations 4.16, 4.32a, and 4.32b, the principal stress state is easily transformed from (σ1 , σ2 , σ3 ) to (ηH , LP , σ ¯H ). The following relation can be derived from Equation 4.32a: 1 σ1 − σ3 =p 2 2 2 σ ¯H F f + Gg + Hh + 2 (Ll2 + M m2 + N n2 )

(4.34)

From the definition of LP and ηH in Equations 4.16, 4.32a, and 4.34, the following relation is obtained:

29

LP σ1 + σ3 = 2ηH − p 2 2 2 σ ¯H 3 F f + Gg + Hh + 2 (Ll2 + M m2 + N n2 )

(4.35)

Using Equations 4.34 and 4.35, the principal stress state is finally expressed in the space of (ηH , LP , σ ¯H ) as follows:   3 − LP σ ¯H σ1 = σm + s1 = ηH + √ 6 T " # 3 − LP = ηH + p σ ¯H 6 F f 2 + Gg 2 + Hh2 + 2 (Ll2 + M m2 + N n2 )   2LP σ ¯H σ2 = σm + s2 = ηH + √ 6 T " # 2LP = ηH + p σ ¯H 6 F f 2 + Gg 2 + Hh2 + 2 (Ll2 + M m2 + N n2 )   3 + LP σ3 = σm + s3 = ηH − √ σ ¯H 6 T " # 3 + LP = ηH − p σ ¯H 6 F f 2 + Gg 2 + Hh2 + 2 (Ll2 + M m2 + N n2 )

(4.36a)

(4.36b)

(4.36c)

where the principal stresses are in the order of σ1 ≥ σ2 ≥ σ3 since the range of the Lode parameter is from −1 to 1. The stress triaxiality defined by the Hill’s 48 equivalent stress can be expressed in terms of (ηv , LP ) which are strictly related to the angles whose combination denotes the direction of principal stress vector in Haigh–Westergaard space as discussed in section 4.2. The angles of (ϕ, θ), which are functions of stress invariants, are graphically represented as shown in Figure 15(a). The stress triaxiality based on the Hill’s 48 equivalent stress has a material dependence as follows: q

L2P + 3 σm σm σ ¯v = = ηv p ηH = σ ¯H σ ¯v σ ¯H 2 F f 2 + Gg 2 + Hh2 + 2 (Ll2 + M m2 + N n2 ) √ √ 2 3I1 cot ϕ = √ ηv = 3 J2

(4.37a) (4.37b)

Consequently, the loading direction is uniquely defined in Haigh–Westergaard space for the anisotropic material by the combination of (ηH , LP ).

4.4

Lode dependence in the principal strain space for anisotropic materials In order to confirm the Lode dependence on the anisotropic material, the vector of the

principal plastic strain rate is derived by the associative flow rule using the Hill’s 48 criterion 30

expressed in terms of the principal stresses. In this derivation, the components of the Rodrigues’ rotation matrix are simply defined as shown in Equation 4.38a for the sake of simplicity in mathematical manipulation:   a b c         R = txy − sz ty 2 + c tyz + sx ≡ d e f      g h i txz + sy tyz − sx tz 2 + c 

tx2 + c

txy + sz txz − sy



    p       ε ˙ A D E σ ε ε ε    1     1  ∂ σ ¯   H = λ˙ =  Dε Bε Fε  ε˙p2 σ2   ˜ ∂ σ          ε˙p  Eε Fε Cε  σ3 3    2  d2 − g 2 Aε           2     Bε  e2 − h2              2  Cε   f 2 − i2 λ˙   =    ¯H  d2 − g 2 e2 − h2   Dε   σ           2      d − g 2 f 2 − i2   E   ε            Fε e2 − h2 f 2 − i2

a2 − g 2

2

    

(4.38b)

   

2 a2 − d2 2d2 g 2 2a2 g 2  2 b2 − e2 2e2 h2 2b2 h2 2 2 2 c −f 2f 2 i2 2c2 i2   a2 − d2 b2 − e2 2degh 2abgh   a2 − d2 c2 − f 2 2df gi 2acgi   b2 − e2 c2 − f 2 2ef hi 2bchi

2 b2 − h2  2 c2 − i2   a2 − g 2 b2 − h2   a2 − g 2 c2 − i2   b2 − h2 c2 − i2

(4.38a)

         2b2 e2       2 2 2c f      2abde        2acdf        2bcef 2a2 d2

 F        G      H   L       M       N

(4.38c)

Substituting Equation 4.16 into Equation 4.38b with the postulate of plastic incompressibility, the following relationship can be described: ε˙p1 Aε +

Aε +

Dε (LP + 1) 2

ε˙p2

= Dε +

Bε (LP + 1) 2

ε˙p3

= Eε +

Fε (LP + 1) 2

Dε Bε Fε (LP + 1) : Dε + (LP + 1) : Eε + (LP + 1) = Pε : Qε : Rε 2 2 2

(4.39a)

(4.39b)

where the principal plastic strain rates are in the order of ε˙p1 ≥ ε˙p2 ≥ ε˙p3 . Besides, the Lode parameter in Equation 4.22 can be also expressed in terms of the principal plastic strains for the anisotropic material as follows: LP =

2Eε ε˙p2 − Fε ε˙p1 − Dε ε˙p3 Fε ε˙p1 − Dε ε˙p3

(4.40)

In the case that the vector of a rotation axis coincides with z-axis in the reference coordinate system, the Lode parameter defined in Equation 4.40 can be written as the following form: 2

2 F sin θS + LP =

 (2N − F − G) 2 2 sin 2θS + Hcos 2θS ε˙3 S ε˙2 − S + Gsin θS ε˙1 − 4    (2N − F − G) 2 2 2 2 F cos θS + Gsin θS ε˙1 − sin 2θS + Hcos θS ε˙3 4

Gcos2 θ



F cos2 θ

2



31



(4.41)

or, equivalently,     r45 (r0 + r90 ) 2 2 r0 sin2 θS + r90 cos2 θS ε˙2 − r0 cos2 θS + r90 sin2 θS ε˙1 − sin 2θS + r0 r90 cos2 2θS ε˙3 2   LP =  r (r + r ) 45 0 90 2 2 r0 cos2 θS + r90 sin θS ε˙1 − sin 2θS + r0 r90 cos2 θS ε˙3 2

(4.42)

The Hill’s 48 equivalent plastic potential, which is the corresponding work-conjugate of the Hill’s 48 equivalent stress, is expressed in terms of the plastic strain rates as s ε¯˙pH =

2

2

2

2

2

2

2(ε˙p23 ) 2(ε˙p13 ) 2(ε˙p12 ) F (ε˙p11 ) + G(ε˙p22 ) + H(ε˙p33 ) + + + F G + GH + HF L M N

(4.43)

With the assumption that there is no shear strain, Equation 4.43 reduces to s ε¯˙pH =

2

2

F (ε˙p1 ) + G(ε˙p2 ) + H(ε˙p3 ) F G + GH + HF

2

(4.44)

Based on the definition of the Hill’s 48 equivalent plastic strain rate as in Equation 4.44, the principal plastic strain rates are finally defined as below:     p        s ε ˙ P   1     ε   F G + GH + HF p ε¯˙ = ε˙p2 Qε     F Pε2 + GQ2ε + HRε2 H       ε˙p     R   ε 3

(4.45)

From the theoretical derivation above, it can be clearly concluded that the principal plastic strain rates are solely dependent on the Lode parameter and the material anisotropy.

32

Chapter 5 Development of an Anisotropic Ductile Fracture Criterion 5.1

Construction of an uncoupled anisotropic ductile fracture criterion As one of the uncoupled approaches to describe the influence of material anisotropy on the

fracture strain, a linear transformation can be considered in modeling the anisotropic ductile fracture criterion, which is discussed by Luo et al. [123] and Gu and Mohr [77] in detail. In this approach, the components of the linear transformation matrix are regarded as the parameters to be calibrated for the fracture criterion. The employment of the modeling approach using the transformation matrix, therefore, can give an increase in the model flexibility in evaluating the fracture strain at a certain loading state considering the material anisotropy. However, the enhancement of the model performance is not always guaranteed when employing the linear transformation in the model construction because it also forces to increase in the non-linearity between the stress state and its corresponding fracture strain predicted from the fracture criterion developed. Each component in the transformation matrix plays a role as a factor to give more weight in certain loading states considering the influence of anisotropy on the fracture initiation. One of the limitations of this approach is that each weight factor is strongly related to the strain or stress components: therefore, the weight factors are determined by the least square scheme which will eventually result in the anisotropic fracture prediction in an average sense. To release the dependence of weight term on the stress or strain components, the components of the linear transformation matrix can be set as the value dependent on the loading direction. There remains, however, a difficulty to make the relationship of each component according to the loading direction in addition to a problem of the high non-linearity between the model parameters. Therefore, there arises a strong need to develop a new approach to deal with the influence of material anisotropy on the fracture initiation. In general, the ductile fracture in metals and alloys is described by the growth and coalescence of voids nucleated at micron scale particles with a large amount of plastic deformation, which eventually causes the loss of the load carrying capacity. Although the successive reduction in load carrying capacity eventually brings about the material softening, it is sometimes useful to develop the ductile fracture criterion as a damage variables which grows during the 33

plastic deformation. In this approach, the damage variable does not affect the constitutive law and the fracture will initiate when the damage variable reaches a certain critical value. This so-called uncoupled ductile damage criterion has a basic form of Z

ε¯pf

0

Z p    1 ε¯f  ˜ ε, ˜ε ˜˙ , . . . d¯ ˜ ε, ˜ε ˜˙ , . . . d¯ εp ≤ C → D = εp ≤ 1 f σ, f σ, C 0

(5.1)

where ε¯pf denotes the equivalent plastic strain at the onset of fracture and C is a material parameter. The damage D can be thought of as dependent on internal variables such as stress, strain, strain rate, temperature, and so on. The Lou–Huh ductile fracture criterion is suggested with the observation of the micro-mechanism for the ductile fracture and consists of three parts denoting damages predicted in a material during deformation. The Lou–Huh ductile fracture criterion has a form of 

2τmax σ ¯v

C1 

1 + 3ηv 2

C2

ε¯pf = C3

 x when x ≥ 0 where hxi =  0 when x < 0

(5.2)

for proportional loading cases [118, 119]. When the loading mechanism is not proportional, the above criterion needs to be evaluated according to the variation of the path with an integral form of Z 0

ε¯pf



2τmax σ ¯v

C1 

1 + 3ηv 2

C2

d¯ εp = C 3

(5.3)

The void nucleation is considered as a function of the equivalent plastic strain, the void growth is described by the stress triaxiality, and the void coalescence is controlled by the normalized maximum shear stress. The governing factors mentioned above are strongly coupled together in the nucleation, the growth, and the coalescence of voids during deformation in the real physical phenomenon. The Lou–Huh ductile fracture criterion can be transformed into the space of the stress triaxiality ηv , the Lode parameter LP , and the equivalent plastic strain ε¯pf at the onset of fracture, as defined in Equation 5.4: 2 q L2P + 3

!C1 

1 + 3ηv 2

C2

ε¯pf = C3

(5.4)

The general stress states can be represented by the combination of ηv and LP , such that the Lou–Huh ductile fracture criterion is capable of predicting the equivalent plastic strain at the onset of fracture over a wide range of stress states. The values of the stress triaxiality ηv and the Lode parameter LP denoting the typical stress states are summarized in Table 13 for isotropic case. 34

For consideration of the anisotropy effect on the fracture strain, the Hill’s 48 criterion is used in combination with the stress triaxiality proposed in section 4.3 to modify the Lou–Huh ductile fracture criterion as follows: 

2τmax σ ¯H

C1 

1 + 3ηH 2

C2

ε¯pf = C3

(5.5)

The anisotropic stress triaxiality proposed is expressed in terms of the principal stresses and the maximum principal stress direction θS through the tensor transformation. With this indicator, the influence of directionality in typical stress states can be considered on the different orientation of sheet metal. The new combination of ηH and LP is hence varied according to the maximum principal stress direction in order to denote the typical stress states as shown in Figure 18(a). In this graphical representation, there exist three stationary states of the uniaxial compression, pure shear, and equi-biaxial tension for the values of ηH , which are distinguished by solid circular symbols as shown in Figures 18(a), 18(b), and 18(c). These values for ηH are theoretically determined under the plane stress condition in the principal stress states and are independent of the maximum principal stress direction. When the material anisotropy is relatively large, the change of ηH also becomes large according to the maximum principal stress direction as shown in Figures 18(b) and 18(c). Here, the r-values of a conventional cold-rolled steel sheet of SPCC 1.0t obtained by Huh et al. [93] are utilized to describe the value of ηH in the case that the material has a large anisotropy. The ηH and LP can be also expressed as functions of the stress path, respectively, under the plane stress condition as illustrated by Figures 18(c) and 18(d). From this theoretical conversion, it is able to calibrate the parameters of the anisotropic ductile fracture criterion directly using the experimental data based on the relation between the stress and strain paths. The change of the value of ηH influences on the construction of the anisotropic ductile fracture criterion and its damage accumulation. The new combination of ηH and LP is summarized in Table 14. The proposed fracture criterion can be represented in the space of the stress triaxiality defined by the Hill’s 48 criterion and the Lode parameter as shown in Equation 5.6, which is derived by substituting Equations 4.16 and 4.32 into Equation 5.5. As a result, the dependence of the equivalent plastic strain to the fracture on the ηH and the LP can be easily evaluated over a wide range of stress states for the anisotropic material. √ −C1  1 + 3η C2 H T ε¯pf = C3 2

35

(5.6)

The anisotropic ductile fracture criterion proposed above can reduce to the Lou–Huh ductile fracture criterion in the case that a material is isotropic. In order to consider a cut-off value for the stress triaxiality depending on the Lode parameter and materials, Lou et al. [120] suggested the Lou–Huh ductile fracture criterion to have a form of 2

!C1

q L2P + 3

  

*

 3 − LP + C2 ηv + q + Co  p 3 L2P + 3  ε¯f = C3

(5.7)

1 + Co

where Co represents the sensitivity of the cut-off value for the stress triaxiality to microscopic structures and temperature because material characteristics are influenced by the change of those factors. In the same manner, the Lode dependence of the cut-off value can be utilized in the anisotropic ductile fracture criterion as shown in Equation 5.8 and the modified cut-off region is represented in Figure 19. * +C2 3 − LP √ −C1 √ + C η + o  H  T ε¯pf = C3 6 T 1 + Co

(5.8)

There are various approaches to determine the value of Co . In the case that it plays an important role in engineering application of the material, the cut-off value can be determined by uniaxial tension test of round bars under the hydrostatic pressure or it can be set to a certain value such that all stress triaxiality at fracture tests are higher than the cut-off value predicted by Equation 5.9. The predicted cut-off value is determined when the numerator of the second term on the left-hand side in Equation 5.8 is equal to zero as follows: ηH +

5.2

3 − LP √ + Co = 0 6 T

(5.9)

Determination of the parameters of the proposed fracture criterion For calibration of the parameters in the anisotropic ductile fracture criterion with a change-

able cut-off value for the stress triaxiality, the constant offset of CO is simply fixed to 1/3 as suggested by Lou et al. [120] since all stress triaxiality at all fracture tests conducted need to be higher than the cut-off values for the stress triaxiality. In order to construct the anisotropic ductile fracture criterion, the three model parameters need to be determined from at least three different equivalent plastic strains at the onset of fracture corresponding to their loading states. 36

The two-dimensional Digital Image Correlation (DIC) method is adopted to measure the strain histories on the surface of three different types of specimens during deformation to calibrate the model parameters directly using the test results. However, the designated stress state induced by the shape of a specimen may not be guaranteed for its own loading state to be perfectly maintained once the material undergoes the severe plastic deformation during the test. For that reason, the effect of the loading path change has to be taken into consideration to evaluate the material behavior properly when the model parameters are calibrated. Recently, Cortese et al. [53] introduced a general damage accumulation law so as to consider the non-proportional loading condition. In this damage accumulation law, the fracture prediction is fulfilled by introducing a damage indicator D as a quantitative measure for ductility consumed at a material point. The range of the damage indicator runs from 0 to 1 and paths of damage growths can vary according to loading conditions and material properties. The damage accumulation law proposed by Cortese et al. [53] has an integral form of

Z D= 0

ε¯p

m  q+1 p ε¯f

m !  p q −1 ε¯f ε¯p d¯ εp ≤ 1 p ε¯f

(5.10)

where q and m are the material parameters which characterize the non-linearity in the damage accumulation and ε¯pf represents the equivalent plastic strain at the onset of fracture predicted from the fracture criterion. This non-linear damage accumulation law coincides with the one proposed by Xue [182] and Papasidero et al. [140] when q = 0. If m = 1 and q = 0, this approach reduces to the linear damage accumulation law. In aid of this damage accumulation law, we can evaluate the onset of fracture under non-proportional loading by performing a structural analysis via numerical simulation. In the present study, the linear damage accumulation law is employed in the process of calibrating the model parameters for simplicity. Under the linear damage accumulation law, Equation 5.8 is expressed as 1 D (¯ ε )= C3 p

Z 0

ε¯pf

C2 n √ −C1  3  X 3 − LP 1 p d¯ ε ≈ di T ηH + √ + 4 3 6 T i=1

(5.11)

Here, the damage indicator D can be evaluated by the sum of incremental damage index di as defined in Equation 5.12a. 1 p −C1 di = Ti C3

  C2 3 3 − (LP )i 1 √ (ηH )i + + [(¯ εpH )i − (¯ εpH )i−1 ] 4 3 6 Ti

37

(5.12a)

where 2 2  1 + (LP )i 1 + (LP )i 2 2 2 2 cos (θS )i + sin (θS )i + G cos (θS )i + sin (θS )i Ti = F 2 2 (5.12b)    N 1 − (LP )i 2 2 2 + H cos 2(θS )i + sin 2(θS )i 2 2 

In Equation 5.11, n denotes the total number of the data for strain path obtained from the test. In order to calibrate the parameters of the above fracture criterion, the principal stress state at the material point where the fracture initiates should be identified. It, however, cannot be directly obtained from the test whose result is the history of the strain field on the surface of a specimen. For the utilization of the strain history on the surface of the specimen obtained from the test as input data to calibrate the parameters of the fracture criterion, the amount of necking needs to be negligible. In the tensile test, however, the triaxial stress state may severely develop in the specimen resulting in the increase of the difference in the level of the localized equivalent plastic strain to the fracture predicted between that in the specimen and that on the surface of the specimen. This will lead to a conservative prediction on the fracture initiation. In the case of AHSS sheets, little amount of necking and thickness reduction is observed before the fracture initiation due to their poor ductility, such that no significant triaxial stress state arises in the specimen when compared with that of mild steels. This can lead to an assumption that there is a slight difference in the level of localized equivalent plastic strain to the fracture predicted between that in the specimen and that on the surface of the specimen. From this assumption, the material point on the surface of the specimen may be expected as an initiation point of the onset of fracture, and thus can be regarded as an alternative one to obtain the data for the construction of the fracture criterion. The data for strain path obtained from the test are the history of strain path defined as a ratio of the minor to major strains under the plane stress condition. The relation between the strain and stress paths, therefore, needs to be defined to calculate the ith component in Equation 5.12b because the stress triaxiality and the Lode parameter can be formulated as functions of the stress path under the plane stress condition. The relation between the strain and stress paths can be derived from the associative flow rule with the Hill’s 48 criterion through β=

A − αB αA − C

38

(5.13a)

where 1 A = (F + G) sin2 θS − H(cos 4θS + 1) − N sin2 θS 2 B = 2(G cos4 θS + F sin4 θS ) + H(cos 4θS + 1) + N sin2 θS

(5.13b)

C = 2(F cos4 θS + G sin4 θS ) + H(cos 4θS + 1) + N sin2 θS α=

ε˙2 ε˙1

(5.13c)

Here, β is a parameter denoting a ratio of the minor to the major stresses under the plane stress condition and α is a parameter denoting a ratio of the minor to major principal strain rates. The anisotropic ductile fracture criterion now can be directly constructed using the history of strain path obtained from tests even though the stress triaxiality and the Lode parameter are defined based on the stress state, which reveals that the stress state should be known to calibrate the original fracture criterion to construct the strain-based fracture criterion. The parameters of the anisotropic ductile fracture criterion are calibrated by minimizing the sum of deviations. The deviation for each loading case is calculated by subtracting the accumulated damage from the unity. Representative strain path and the calibrated model parameters are utilized for calculation of the damage accumulation through the fracture criterion. The algorithm to optimize the parameters of the fracture criterion is summarized as follows: min Ci

M X N X k 1 − D (η , L , θ ) P S j H

(5.14)

k=1 j=1

The subscript j denotes the typical loading cases determined for the fracture tests such as pure shear, uniaxial tension, and plane strain tension test and N denotes the total number of typical loading cases. The superscript k and M denote the orientation of sheet metal and the total number of the orientations determined for the fracture tests, respectively. For instance, when k = 1, it represents the difference evaluated from the test results for the rolling direction °

0°, that is, δj1 ≡ δj0 . In order to minimize the above objective function, the Nelder–Mead simplex algorithm built in the MATLAB program is used instead of the iterative algorithms. In general, the iterative algorithms have solid mathematical foundations; however, the complicated derivation of iterative algorithms makes the application of those in many problems being difficult. Even if the derivation can be obtained, the results from the optimization strongly depend on the initial values and an incorrect choice of initial values may lead to severe local minima. Relatively, the heuristic algorithms depend less on the properties of the objective functions and can provide approximate solutions with robust stabilization despite the lack of strict mathematical descriptions and expensive computational cost. Due to the high nonlinearity of the 39

objective function, the Nelder–Mead algorithm is suitable for dealing with the current problem. The obtained parameters of the proposed fracture criterion are summarized in Table 15. The three-dimensional fracture envelopes predicted for the three loading directions are represented in Figure 20. A fracture locus is also extracted from the three-dimensional fracture envelope by projecting the fracture strains predicted under the plane stress condition onto the space of the stress triaxiality for application to sheet metal forming processes. From the fracture loci constructed at the three loading directions of 0°(RD), 45°(DD), and 90°(TD), the anisotropic ductile fracture criterion proposed can take account of the equi-biaxial fracture strain which is independent of the loading direction. It is simply confirmed that the equi-biaxial fracture strain is kept a constant value regardless of the variation of the maximum principal stress direction as shown in Figure 21. The damage evolutions are computed for all loading paths tested at the three loading directions as shown in Figure 22 to evaluate the equivalent plastic strain to the fracture predicted from the anisotropic ductile fracture criterion constructed. Through this computation, it is evaluated that equivalent plastic strains to the fracture are predicted well for the rolling direction meanwhile those are somewhat overestimated at the loading conditions of uniaxial tension and plane strain tension for the diagonal and transverse directions. These results may relate to an optimization algorithm employed and the number of loading paths utilized for calibration of the model parameters of the anisotropic ductile fracture criterion. It is, consequently, necessary to verify various optimization methodology and what kinds of loading paths should be utilized for the model construction to properly evaluate the equivalent plastic strain to the fracture over a wide range of loading paths. There is one purely empirical way to accomplish the better prediction on the equivalent plastic strain to the fracture and it will be discussed in section 7. As a result of the construction of the anisotropic ductile fracture criterion, damages of each element in FE simulation can be properly evaluated in consideration of the material anisotropy which may vary with the plastic deformation. This approach, however, needs so tremendous computing time to predict the fracture initiation in numerical simulation of complicated forming processes. It is, therefore, hard to be applied to efficiently and effectively predict the formability of the material. As a compromise, a representative fracture-based forming limit criterion is proposed as a reference to predict the onset of fracture for the anisotropic material in complicated sheet metal forming processes.

40

Chapter 6 Construction of Fracture Forming Limit Diagrams 6.1

Strain-based fracture forming limit diagram The strain-based Forming Limit Diagram (FLD) is conventionally utilized in the sheet metal

forming community to predict the forming severity. It is expressed in the space of the major and minor strains by assuming proportional loading states under the plane stress condition. In sheet metal forming processes, the conventional FLD with the safety margin can sufficiently guarantee successful products without a defect such as a localized thinning which leads to a tear in sheet metals. Prediction of the onset of fracture in AHSS sheets needs a Fracture Forming Limit Diagram (FFLD) to evaluate the forming limit in sheet metal forming processes since AHSS sheets fail with a very little amount of necking as discussed by Li et al. [113] and Beese et al. [27]. The FFLD as an alternative to the necking-based FLD is derived by equating the fracture strain obtained from the fracture criterion with the equivalent plastic strain corresponding to the applied yield criterion. The equivalent plastic strain of the Hill’s 48 criterion can be expressed as Equation 6.1 which is the corresponding work-conjugate of the Hill’s 48 equivalent stress. ε¯˙p =

ε˙p1

r

F + Gα2 + H(1 + α)2 F G + GH + HF

(6.1)

The range of α is from −2 to 1, which denotes the strain states from the uniaxial compression to the equi-biaxial tension. The major and minor strains at the onset of fracture now can be evaluated from Equation 6.2 under the proportional loading condition. εp1

=

ε¯pf



F + Gα2 + H(1 + α)2 F G + GH + HF

− 12

, εp2 = αεp1

(6.2)

where ε¯pf denotes the fracture strain predicted from the strain-based fracture criterion. Figure 23 represents the FFLDs obtained according to the three loading directions. In these representations, typical strain paths are varied according to the tensile loading direction which is the angle inclined from the rolling direction of sheet metal. These results are obtained from the variation of the maximum principal stress direction with the material anisotropy and it can explain the anisotropic deformation behavior at different orientations of sheet metal under 41

the typical stress states. For the equi-biaxial tension condition, the major and minor strains predicted at the onset of fracture are varied according to the loading direction. This explains that the fracture can actually initiate along any direction in sheet metal although the equivalent plastic strain to the fracture is independent of the loading direction. The FFLDs constructed demonstrate that there is no general tendency for the predicted fracture strains of typical loading paths for the three loading directions because the data utilized for the calibration of the model parameters show a complicated tendency on the equivalent plastic strain to fracture. FFLDs predicted from different directions of the maximum principal stress are also reviewed for comparison to evaluate the formability of anisotropic materials as shown in Figure 24. This provides direct information for predicting the overall formability of AHSS sheets with the material anisotropy in a design process for manufacturing. In comparison of the fracture-based forming limits suggested, it can be recommended that the material formability should be considered by the lowest forming limits over a wide range of loading conditions because the direction of fracture initiation cannot be predicted in complicated forming processes. The strain-based FFLD is theoretically constructed under the proportional loading condition. The fracture limits, however, might be expected to have a path dependence for other loading paths and more complicated forming processes. As one of the approaches to deal with this issue, the concept of path-independent fracture forming limits is suggested in section 6.4.

6.2

Scaling of the strain-based fracture forming limit diagram according to the direction of major stress axis As discussed in the previous section, it is necessary to construct a representative fracture-

based forming limit curve for the anisotropic material to enhance computing time efficiency and capability of the evaluation for the onset of fracture before manufacturing of a product designed. Stoughton et al. [163, 166] discussed the construction of a single forming limit to monitor the failure for anisotropic materials. In this section, a solution to the challenge of constructing a single fracture-based forming limit curve for the anisotropic material is suggested based on the concept of FLD for anisotropic materials discussed by Stoughton et al. [163, 166]. A strain-based FFLD for the rolling direction is regarded as a reference forming limit to scale the FFLD with various orientations of a sheet metal. The reference forming limit can be defined as      − 12   n o εF  2 2 1 F + Gα + H(1 + α) 1 = ε¯pf εF = εF  α F G + GH + HF 1 42

(6.3)

Based on the concept of a generic anisotropic forming limit discussed by Stoughton et al. [163, 166], the major and minor strains at the onset of fracture can be suggested as follows:       1 1 n o εF  1 = εF1 (α, θS ) = εF1 (α)fI (θS ) (6.4a) εF = α α εF  1

εF1 (α, θS ) =

  εF1 (α)fI (θS )

α ≤ αpl 

α − αpl αeb − α feb  fI (θS ) + εF1 (α)fII (θS ) = εF1 (α) RD αeb − αpl αeb − αpl feb

 α > αpl

(6.4b)

with fI (θS ) =

i 1 h 2 4 4 2 f cos θS + (4f45 − f90 − f0 ) cos θS sin θS + f90 sin θS f0 0 αeb =

F F −H , αpl = G 2G + H

(6.4c) (6.4d)

where fI (θS ) is the normalized fracture forming limit in the plane strain and θS is the angle of the major stress to the rolling direction of sheet metal, which can be also interpreted as the maximum principal stress direction. The values of (f0 , f45 , f90 ) are the FFLD0 predicted from the anisotropic ductile fracture criterion along the rolling, the diagonal, and the transverse RD indicate the forming limit under the equi-biaxial tension and the equidirections. feb and feb

biaxial forming limit predicted from the FFLD for the rolling direction, respectively. Typical strain paths can be theoretically calculated based on the associative flow rule under the plane stress condition. Here, (αeb , αpl ) indicate the theoretical strain path in the equi-biaxial tension and that in the plane strain based on the Hill’s 48 criterion, respectively. Since the forming limits predicted under the equi-biaxial tension should be the same, the right-hand side of FFLD is scaled using a scaling factor, fII . Figure 25 shows the FFLD0 according to the maximum principal stress direction and Figure 26 represents the scaling factor in the space of the maximum principal stress direction and the strain path. With the employment of the approach suggested in Equation 6.4a, the FFLDs scaled according to the maximum principal stress direction are constructed and compared with the FFLDs predicted from the anisotropic ductile fracture criterion as represented in Figure 27. The FFLD for the rolling direction is evaluated as the lowest forming limit among FFLDs scaled because the level of FFLD scaled is determined based on the level of FFLD0 . The maximum principal stress direction for the lowest FFLD scaled would be varied if the combination of values of FFLD0 changes influencing on the curve for FFLD0 predicted. In comparison of the FFLDs predicted from the anisotropic ductile fracture criterion and the FFLDs scaled, the FFLD for the rolling direction is located at the lowest region showing that it can be considered as a representative fracture-based forming limit. 43

6.3

Stress-based fracture forming limit diagram The conventional FLD constructed under the proportional loading condition has a limi-

tation that it can be only utilized for the situation when the strain path is linear or nearly linear during forming processes to make a reliable prediction for forming severity. Therefore, the safety margin is usually taken into consideration for the complex forming process where the strain path is non-linear. In contrast, the stress-based FLD has little dependence on the strain path change. It means that when the strain-based FLD is transformed into the stress-based FLD, the intrinsic dynamic nature of the strain-based FLD degenerated in the stress space as discussed by Arrieux et al. [5], Stoughton [161], and Stoughton and Yoon [163,166]. The stressbased FFLD can be constructed from the stress-based three-dimensional fracture envelope in a similar way to the construction of the strain-based FFLD. The equivalent stress at the onset of fracture can be obtained based on an adequate hardening model which accurately describes the material behavior. The utilization of the adequate hardening model plays a great role in the computation of the stress states for damage calculation because the hardening model directly affects the stress states at each element in FE analysis. In order to precisely predict the stress state in finite element analysis, the Gurson–Tvegaard–Needleman (GTN) material model could be utilized to predict the stress states after necking, which is beyond of the scope of the present study. In this study, the hardening curves obtained from experiments is utilized with the mixed Swift–Voce hardening model given in Equation 3.1 to compute the stress states appropriately although the range of tolerance should be taken into consideration for the stress states after necking. For the theoretical transformation from the strain-based FFLD to the stress-based FFLD, the critical fracture stress is predicted by the mixed Swift–Voce hardening model in combination with an anisotropic hardening model suggested by Stoughton and Yoon [165] to include anisotropic behavior into the stress-based FFLD, i.e.,



σ ¯θ θS , ε¯pf



v −1  u u cos2 θ 2 2 sin θS  sin 2θS u S   cos 2θS +   = t   − p p 2 2 σ ¯02 ε¯f σ ¯90 ε¯f σ ¯45 ε¯pf

(6.5)

where σ ¯0 , σ ¯45 , and σ ¯90 are the mixed Swift–Voce hardening curves fitted at the three loading directions. Here, θS and ε¯pf denote the maximum principal stress direction and the equivalent plastic strain at the onset of fracture predicted from the fracture criterion, respectively. Following this methodology, the stress-based three-dimensional fracture envelopes were constructed and represented in Figure 28. The stress-based FFLD is theoretically constructed by imposing the equivalence between the critical fracture stress predicted by Equation 6.5 and the Hill’s 48 44

equivalent stress under the plane stress condition as follows:     1 n o σ F  p −1 1 F = =σ ¯θ Tβ σ β  σ F  1

(6.6a)

where 2 2  1 + LP (β) 1 + LP (β) 2 2 2 2 Tβ = F cos θS + sin θS + G cos θS + sin θS 2 2    1 − LP (β) 2 N 2 2 + H cos 2θS + sin 2θS 2 2 

(6.6b)

Here, LP can be expressed in terms of the stress path β according to typical stress states as shown in Figure 18(d). The explicit relations between LP and β are defined as below: σ2 σ1 σ3 β= σ1 σ1 β= σ3 σ2 β= σ3

  2β − 1      −1 − β    1−β LP (β) = −1 −β     β−1     −2β + 1

β=

∈ [0, 1]

(σ1 > σ2 > 0, σ3 = 0)

∈ [−1, 0] (σ1 > 0, σ2 = 0, σ3 < 0, |σ1 | > |σ3 |) (6.6c) ∈ [−1, 0] (σ1 > 0, σ2 = 0, σ3 < 0, |σ3 | > |σ1 |) ∈ [0, 1]

(σ1 = 0, 0 > σ2 > σ3 )

The stress-based FFLDs obtained according to the three loading directions and the corresponding fracture loci are represented in Figure 29 dealing with the stress states ranging from the equi-biaxial compression to the equi-biaxial tension for application to complicated sheet metal forming of AHSS sheets. The stress-based FFLDs predicted from various directions of the maximum principal stress are represented in Figure 30(a) and the stress-based FFLD corresponding to the loading direction of 45° is evaluated as the lowest forming limit in this comparison. Due to the characteristic of the anisotropic hardening model given in Equation 6.5, the equi-biaxial fracture stresses, however, solely depend on the loading direction, which is not plausible in light of the non-directionality of the equi-biaxial stress state. One point convergence should be preserved at the equi-biaxial stress state. For the one point convergence of the equi-biaxial fracture stress regardless of the variation of loading direction, the fracture stresses are then normalized by the equi-biaxial fracture stress at each loading direction:     n o σ F  p −1  1  σ ¯θ 1 Tβ = f σF = f β  σ F  σ ¯ (θ , ε ¯ ) S b b 1

(6.7a)

where 

σ ¯fb θS , ε¯bf



v −1  u u cos2 θ 2 2 sin θS  sin 2θS u S   cos 2θS +   = t   − 2 2 σ ¯02 ε¯bf σ ¯90 ε¯bf σ ¯45 ε¯bf 45

(6.7b)

Here, σ ¯fb and ε¯bf denote the equi-biaxial fractre stress and strain. Figure 31 represents the stressbased FFLD in consideration of the one point convergence at the equi-biaxial stress state.

6.4

Polar effective plastic strain-based fracture forming limit diagram The stress-based FLD has a demerit that it is difficult to quantify the margin of safety

even though it has a merit of path independent nature. Significant changes in strain take place at stress levels close to the necking limit due to the large strain increment with small stress increment in the true stress–strain curve, such that it forces to increase an effort on accurate measurement of the strain in the necking region. It is, however, not easy to measure the strain in the necking region and there still arises a question of reproducibility in measuring the strain even though the strain evaluation is conducted in the necking region. To overcome the drawback of the stress-based FLD, the forming severity is evaluated based on the equivalent plastic strain as discussed by Zeng et al. [191]. The path dependence of the strain-based FLD degenerated into one single curve when it is converted into the space of the equivalent plastic strain with respect to the stress path or the strain path. These approaches have advantages that it is relatively easy to confirm the margin of safety near the necking limit and there is no dependence on the stress–strain relation, such that the FLD might be extendable to various material models which have a non-monotonic relation between the stress and strain. Recently, Stoughton et al. [167] have proposed a mathematically equivalent solution to using the variables of the equivalent plastic strain and the strain path as proposed by Zeng et al. [191] in a polar diagram, which might be more appealing for application to automotive industries than the previous FLDs. The conventional FLD in polar coordinates, which is called Polar Effective Plastic Strain (PEPS) Forming Limit Diagram (FLD), has a path-independent nature which is supported not only by experimental evidence but also by theoretical analysis in the works by Stoughton and Zhu [168] based on bifurcation analysis, and Zeng et al. [191] based on analysis with Marciniak–Kuczynski (M–K) model. In this study, it is assumed that a polar effective plastic strain (PEPS) fracture forming limit diagram (FFLD) has a path-independent nature when the strain-based FFLD is converted into a polar diagram which is the transformation of the path dependent diagram to the path independent space similar to the path independence of the PEPS FLD. The conversion of the FLD into the polar diagram utilizes each strain path on the forming limit curve and corresponding equivalent plastic strains because the strain path can be defined

46

as the arctangent of the ratio of the principal strain increments as follows: 0

−1

θ = tan



∆εp1 ∆εp2

 (6.8)

The variables of the PEPS FLD in a Cartesian equivalent system can be then obtained as (x, y) = (¯ εp sin θ0 , ε¯p cos θ0 )

(6.9)

where ε¯p is the radial variable in the polar diagram and the direction θ0 reflects the strain path in the Cartesian coordinate system superimposed on the polar diagram. By applying this procedure to the strain-based FFLD constructed under the proportional loading condition, it can be converted into the PEPS FFLD which is assumed to be path independent. In this representation, the variation of loading condition prior to reaching the critical value is considered by plotting the deformation history of a material in the Cartesian equivalent system with the PEPS FFLD together as considered in the PEPS FLD. The PEPS FFLDs obtained are represented in Figure 32 and it is clearly observed that the onset of fracture with the PEPS FFLDs shows notable difference among the three loading directions especially in the shear-dominated region. It reveals that the anisotropy strongly affects the equivalent plastic strain to the fracture. Consequently, the fracture prediction has to be dealt with the lowest equivalent plastic strains to the fracture predicted from the fracture criterion with the consideration of the maximum principal stress direction. In comparison of the PEPS FFLDs converted from the strain-based FFLDs and the ones from the FFLDs scaled, the PEPS FFLD converted from the strain-based FFLD for the rolling direction is evaluated as the lowest forming limit as shown in Figure 33. Hence, the PEPS FFLD for the rolling direction can be regarded as a representative fracture forming limit curve for the DP980 1.2t steel sheet to capture the onset of fracture in complicated forming processes in the case that the non-proportional loading condition is considered for the calibration of the parameters of the proposed fracture criterion.

47

Chapter 7 Emprical Way of Modeling the Fracture Criterion Considering the Material Anisotropy 7.1

Model parameters dependent on the maximum principal stress direction In section 5, the anisotropic ductile fracture criterion is suggested by employing the Hill’s

48 criterion on the basis of the Lou–Huh ductile fracture criterion. In the derivation of the anisotropic ductile fracture criterion, the principal stress states are theoretically expressed in the space of (ηH , LP , σ ¯H ). As a result, the proposed fracture criterion has a capability of describing the anisotropy effect on the deformation behavior of sheet metal according to the maximum principal stress direction. There is, however, a little difference in comparison of damage evolution predicted by the proposed fracture criterion as shown in Figure 22 because the number of experimental data utilized in the model construction is greater than that of model parameters. In an attempt to predict the equivalent plastic strain to the fracture more precisely over a wide range of stress states, the model parameters can be suggested as a function of the angle between the maximum principal stress axis and the rolling direction of sheet metal in a purely empirical way, i.e., Ci (θS ) = Pi cos4 θS + (4Qi − Ri − Pi ) cos2 θS sin2 θS + Ri sin4 θS

(7.1)

where Pi , Qi , and Ri denote the model parameters. Note that the directional dependence in Equation 7.1 vanishes when the parameters of Pi , Qi , and Ri are set to the same value. For the equi-biaxial stress state (ηH = 2/3, LP = 1) under the plane stress condition (σ3 = 0), Equation 5.8 with Equation 7.1 gives 

1 √ F +G

C1 (θS )

ε¯bf = C3 (θS )

(7.2)

where ε¯bf represents the equivalent plastic strain at the onset of fracture under the equi-biaxial stress state. It is worth to note that ε¯bf becomes C3 (θS ) only if the equivalent stress is set to the equi-biaxial yield stress, which leads to F + G = 1. Setting the model parameters of P3 , Q3 , and R3 to the value of the equi-biaxial fracture strain subsequently, the non-directionality of the equi-biaxial fracture strain is eventually fulfilled. In a similar analogy to demonstrate the

48

non-directionality of the equi-biaxial fracture strain, a new anisotropic ductile fracture criterion is proposed as * +C2 (θS ) 3 − LP √ −C1 (θS )  ηH + 6√T + CO  T ε¯pf = C3 (θS ) 1 + CO

(7.3)

For the confirmation of enhancement in the fracture predictability when employing the above approach, the equi-biaxial yield stress and the equi-biaxial fracture strain are simply set as the yield stress for the rolling direction and the equi-biaxial fracture strain predicted by the newly developed criterion, respectively. As a consequence, fracture forming limits, as well as the damage accumulation, are re-evaluated as shown in Figures 34 and 35. The obtained model parameters and the anisotropic parameters of the Hill’s 48 criterion are summarized in Tables 16 and 17, respectively. In the comparison of the obtained results with the ones evaluated in sections 5 and 6, it can be simply confirmed that the fracture predictability is significantly enhanced when the above empirical way is employed for the construction of the anisotropic ductile fracture criterion. This enhancement in the fracture predictability is closely associated with an increase in the number of the model parameters as well as an adequate description of the dependence of both the Lode parameter and the anisotropic stress triaxiality. The increase in the number of the model parameters, however, does not always ensure the enhancement of the model performance when the model has a high non-linearity between a set of input variables and its corresponding output values. Meanwhile, the form of the proposed fracture criterion shows a clear linearity when it is reviewed in a logarithmic form as follows:  C1 (θS ) ln

1 √ T



* + C2 (θS ) ln 

+ 3 − LP ηH + √ + CO  + ln ε¯pf = ln C3 (θS ) 6 T 1 + CO

(7.4)

It is thus without a doubt that the enhancement of the model performance can be ensured in the present form of the fracture criterion proposed with the increase in the number of the model parameters.

7.2

Weight function considering non-directionality of the equibiaxial fracture strain As discussed in section 7.1, in the viewpoint of the fracture prediction over a wide range of

stress states, the fracture predictability of the fracture criterion proposed strongly depends on whether the non-directionality of the equi-biaxial strain holds. In general, the parameters of the 49

facture criterion for advanced metals are calibrated by the quantities from typical loading states of the pure shear, the uniaxial tension, and the plane strain tension: those loading states play a significant role in determining the overall shape of the FFLD as well as a three-dimensional fracture envelope. If the non-directionality of the equi-biaxial fracture strain is considered in the calibration of the model parameters by eliminating the directional dependence of the parameter C3 (θS ) through P3 = Q3 = R3 = ε¯bf , there may arise potential loss of the model predictability on the fracture strains for other loading states due to the inflexibility forced to the parameter C3 (θS ) = ε¯bf = C3 . In an effort to overcome this drawback, a weight function is considered with Equation 7.3: * +C2 (θS ) 3 − LP √ −C1 (θS ) √ η + + C H O   w (˜ x) ε¯pf = C3 (θS ) T 6 T 1 + CO

(7.5)

˜ stands for the state variable vector. The weight function w(˜ where x x) is initially defined by the Cauchy stress tensor σ so as to deal with the general stress states. Since the Cauchy stress tensor is directly related to the anisotropic stress triaxiality, the Lode parameter, and the maximum principal stress direction, an equivalent form of the weight function can be expressed as w(ηH , LP , θS ). For the non-directionality of the equi-biaxial fracture strain, a value of the weight function is necessary to become unity under the equi-biaxial stress state. The weight function is thus simply proposed to have a form of !C4 (θS )

1 w (˜ x) ≡ w (ηH ) =

2 + ηH

5 9

(7.6)

With the adoption of the weight function, the final form of the uncoupled anisotropic ductile fracture criterion can be rewritten with the postulate of the plane stress condition for sheet metal forming application as follows: * +C2 (θS ) 3 − LP √ −C1 (θS ) √ η + + C H O   T 6 T 1 + CO

50

!C4 (θS )

1 2 + ηH

5 9

ε¯pf = C3

(7.7)

Chapter 8 Discussion of Anisotropic Fracture Forming Limits for Application to the Automotive Industry In previous chapters, modeling methodologies to predict the fracture limits were rigorously discussed in consideration of the material anisotropy. All kinds of fracture limits were constructed based on the non-proportional loading condition in calibrating the parameters of anisotropic fracture criterion proposed in an integral form. In this section, each fracture forming limit for application to the automotive industries will be re-constructed by the anisotropic fracture criterion discussed in section 7.2. For quantitative evaluation of the fracture predictability over a wide range of the stress states, the strain histories will be excluded such that the fracture strains given in Table 7 will be used as the physical quantities for each loading state designated. This assumption lies in a fact that the strain paths were maintained almost constantly during each test except for the uniaxial tension case: the strain path during the uniaxial tension test slightly changed after the necking as shown in Figure 9. Since the change of strain paths was not severe during the tests, we can simply assume that the material is under the proportional loading condition before the fracture initiates. In this viewpoint, the equivalent plastic strain to fracture evaluated from the test can be regarded as the fracture strain corresponding the targeted loading condition induced by the specimen geometry, which can make possible to quantitatively evaluate the model performance in predicting the onset of fracture by comparing the fracture strains obtained from the tests with the ones predicted from the fracture criterion. Besides, there may arise unreasonable predictions on the fracture strain from the fracture criterion calibrated under the non-proportional loading condition due to possible influences among the strain path change on the damage accumulation which will influence the value of fracture strain predicted at a certain stress state: that is, the change of strain path used in the model calibration will cause an increase in uncertainty of the fracture strains predicted from the fracture criterion for a certain loading state. In this background, the proportional loading assumption can give direct intuition for assessment of the formability of a material in comparison of the fracture strains obtained from the experiments and the ones predicted from the proposed fracture criterion. In the case of the stress-based FFLD, it is essential to obtain the flow stress curve in the

51

large deformation range because the range of fracture strains evaluated from the tests is from about 0.2 to about 1. For this issue, the hydraulic bulge test will be employed to evaluate the stress–strain response in the large deformation range because it is impossible to obtain a reliable relation between the stress and strain from the uniaxial tension test particularly for AHSS sheets. For example, the uniform elongation of the DP980 1.2t in this study was measured as around 0.045.

8.1

Hydraulic bulge tests for quantification of the stress–strain response in the large deformation range The uniaxial tension test is one of the most conventional testing methods to characterize

the hardening behavior of sheet metals. For mild steels, this method can provide the flow stress curve over the plastic strain of about 50%. On the contrary, for AHSS sheets, it is difficult to characterize the flow stress curve in the large deformation range by the uniaxial tension test because AHSS sheets generally exhibit small ductility. In general, the stress– strain response beyond the point of the uniform elongation is predicted by extrapolating the flow stress curve using various mathematical models such as Ludwik, Swift or Gosh, and so on. Due to the characteristics of mathematical models, the flow stress curves extrapolated cannot be represented as a unique curve which provides the inherent material behavior. One, hence, should have a careful consideration in determining the mathematical model to predict the material behavior beyond of the uniform elongation because a trustworthy evaluation of the material behavior can be fulfilled only if the accurate material properties are given in numerical simulation. As one of the ways to obtain the inherent stress–strain response at a higher level of plastic deformation, the hydraulic bulge test incorporating the DIC method has been discussed by many research [104, 115, 130, 131, 134]. The most important advantage of the hydraulic bulge test is the inexistence of frictional interactions between tools and a specimen in the area of interest, which not only simplifies the analytical solutions for the calculation of stress and strain but also ensures the reproducibility of the test.

8.1.1

Theory for stress and strain measurement

The determination of the flow stress curve is based on analysis of measurable parameters obtained from the hydraulic bulge test. Through Laplace’s equation from the membrane theory, it is possible to make a theoretical relationship used for the calculation of the stress at the apex of the dome. The validity of the relationship obtained is mainly related with a ratio of the sheet

52

thickness to the bulge diameter. For a spherical membrane with a very small ratio between the radius of curvature and the polar thickness, the meridian stress is much higher than the bending stress such that the effects of bending can be neglected [6, 149]. Considering an axially symmetric element in a thin sheet, the equilibrium equation under the action of uniform pressure p can be written as σ1 σ2 p + = ρ1 ρ2 t

(8.1)

where σ1 and σ2 are the principal stresses on the sheet surface, ρ1 and ρ2 are the corresponding radii of the curved surface at mid-thickness, and t is the sheet thickness, respectively. For the axisymmetric case of the bulge test, both principal stresses can be regarded as the equivalent one (σ1 = σ2 = σb ) and equal to the membrane stress which is also called biaxial stress. From the same reason, the curvature radii also remain equivalent (ρ1 = ρ2 = ρ) whatever plane is considered. With those conditions, the biaxial flow stress can be determined by σb =

pρ 2t

(8.2)

By evoking the assumption that the material is incompressible and the shape of deformed curvature is spherical, the biaxial strain εb is equal to the true thickness strain εt at the polar region. The true thickness strain is defined as  εb = ln

t t0

 = −(ε1 + ε2 )

(8.3)

where t0 is the initial sheet thickness, ε1 and ε2 are the principal plastic strain components. In the case that the ratio of the sheet thickness to the bulge diameter is not small enough, the bending effect should be taken into account to evaluate the biaxial flow stress curve in more detail [178]:

(1) In the thickness direction, the strains of the outer layer are bigger than those in the middle of the sheet or as the average over the thickness because the strains in the material are influenced by both stretching and bending. This implies that the average strains are overestimated (true strain becomes large) and the actual thickness is underestimated (true stress becomes large). (2) The radius of the surface curvature ρs is bigger than the radius of the middle layer ρm (true stress becomes large). 53

(3) The pressure affects the inner surface of the material only. All theses single influences lead to higher values of both the stress and strain, such that the flow stress curve will be overestimated typically in the standard procedure. This overestimation can be compensated by using the ARAMIS system for the optical measurement. The bending influence in strain determination is eliminated by calculating the strain values and coordinates directly for the middle layer of the sheet. The actual thickness can be calculated by Z t=

1 t0 dt ≈ tm = [1 + e1 (t)] [1 + e2 (t)] [1 + e1 (t)] [1 + e2 (t)]

(8.4)

where e1 (t) and e2 (t) denote the engineering strain distributions in the thickness direction and the subscript m denotes the middle layer of the sheet. Note that the relation in Equation 8.4 is valid only if the strains of e1 (t) and e2 (t) are linearly distributed in the thickness direction and the volume keeps constant during the deformation. The radius of the curvature in the middle layer ρm can be determined by the following equation: ρm = ρos −

tm tm = ρis + 2 2

(8.5)

where ρos and ρis denote the outer and inner radius of the curvature of the sheet as shown in Figure 36, respectively. Equation 8.5 is valid only for the case of a homogeneous thinning in the sheet. If a localized thinning at the apex of the dome takes place, the above relation does not give reasonable results for the real ρm . From a balance of forces caused by the pressure on the inner surface and the tangential stresses in the shell, the equation for biaxial flow stress can be derived as σb (2πρm )tm = pπ ρis

2

(8.6)

Substituting Equation 8.5 into Equation 8.6 gives   tm 2 pρm 1 − 2ρm σb = 2tm

(8.7)

Thus, to characterize the stress–strain response by Equation 8.7, it is needed to obtain two quantities such as the radius of the curvature ρm and the thickness tm which are related to the strains at the apex of the dome as given in Equations 8.4 and 8.5. The engineering strains at the outer layer of the dome, which are used for the determination of ρm and tm can 54

be computed during the bulge test based on the 3D-DIC method in aid of the ARAMIS system. The radius of the curvature at the apex P (x0 , y0 , z(x0 , y0 )) is computed from the 3D-coordinates measured from the ARAMIS system and it is defined as the negative reciprocal value of the mean curvature at the apex as follows [3]: ρos = −

2 2 =− κ1 + κ1 zxx (x0 , y0 ) + zyy (x0 , y0 )

(8.8a)

where the location of the apex of the dome is determined by ∂z (x0 , y0 ) ∂ 2 z (x0 , y0 ) = zx (x0 , y0 ) = 0 and = zxx (x0 , y0 ) < 0 ∂x ∂x2 ∂z (x0 , y0 ) ∂ 2 z (x0 , y0 ) = zy (x0 , y0 ) = 0 and = zyy (x0 , y0 ) < 0 ∂y ∂y 2

(8.8b)

Consequently, the flow stress curve can be obtained on the basis of three variables in aid of the ARAMIS system: the internal pressure recorded during the bulge test; the bulge radius; and the thickness evaluated at the apex of the dome.

8.1.2

Flow stress curves at the three loading directions considering large deformation

The hydraulic bulge test is performed for the DP980 1.2t steel sheet to evaluate the equibiaxial flow stress curve using the testing apparatus in POSCO as shown in Figure 37. Test conditions are briefly summarized in Table 18. The flow stress curves obtained from the hydraulic bulge test are represented in Figure 38. In the experiments, the equi-biaxial fracture strain was evaluated at the material point showing the maximum equivalent plastic strain before the fracture initiation as shown in Figure 39. Since a localized neck initiates around the equivalent plastic strain of 0.45 in the experiments, the flow stress curves were utilized until that of 0.4 for the characterization of the flow stress curve in the large deformation range. Theoretically, the localized neck forms at a characteristic angle at which the incremental strain in that direction, εθ , becomes zero: r 2

2

dεθ = dε1 cos θ + dε2 sin θ = 0 → tan θ =



dε2 dε1

(8.9)

The localized neck cannot be thus generated under the equi-biaxial condition because dε2 is positive, such that there is no angle at which the localized neck forms; however, in practice, a small pre-existing groove perpendicular to the largest principal stress can grow gradually

55

into the localized neck and how rapidly this initiates depends largely on the strain hardening exponent of n and to a lesser extent the strain rate exponent of m. For the characterization of the flow stress curve in large deformation range, the flow stress curve obtained from the hydraulic bulge test was used in combination with that obtained from the tensile test. On the basis of the plastic work equivalence principle as depicted in Figure 40, the flow stress curve after the uniform elongation was predicted by the curve scaled from the flow stress curve evaluated from the hydraulic bulge test. In the way of scaling, the plastic works per volume were first calculated from the flow stress curves of each test as shown in Figure 41(a). After determining the equivalent plastic strain corresponding to a certain value of the plastic work per volume, a ratio of directional yield stress to the equi-biaxial yield stress was then obtained according to the plastic work per volume as shown in Figure 41(b). As a consequence, the flow stress curve obtained from the hydraulic bulge test was scaled by using a normalized yield stress ratio saturated with respect to the plastic work per volume for prediction of the stress–strain response after the uniform elongation. Figure 42 represents the flow stress curves in the large deformation range with the curves fitted by the mixed Swift–Voce hardening model. The model parameters are summarized in Table 19. The mechanical properties obtained from the hydraulic bulge test and the anisotropic parameters of the Hill’s 48 criterion are summarized in Tables 20 and 21. Here, the Hill’s 48 anisotropic parameters are obtained by the saturated values of the normalized yield stress ratios.

8.2

Fracture limits under the proportional loading condition Using the obtained results from various fracture tests including the hydraulic bulge test, the

fracture envelopes were re-evaluated at the three loading directions of 0°, 45°, and 90° under the proportional loading condition as shown in Figure 43. The parameters of the fracture criterion given in Equation 7.7 are calibrated with use of the fracture strains averaged at each loading condition. The identified parameters are summarized in Table 22. As can be seen in Figure 43, the proposed fracture criterion shows not only the characteristic of non-directionality of the equibiaxial fracture strain but also the fracture predictability with great accuracy over a wide range of stress states according to the loading direction. For the confirmation of the model performance on the fracture prediction, fracture loci predicted from the proposed fracture criterion are represented together with those predicted from the anisotropic ductile fracture criterion before adopting the weight function as shown in Figure 44. Note that both criteria employ the model parameters dependent on the maximum principal stress direction which are discussed in section

56

7.1. In this comparison, it is simply confirmed that the fracture predictability is remarkably enhanced especially for the uniaxial tension condition including various loading states when the weight function is employed. This enhancement in the model prediction can be also confirmed from the strain-based FFLDs as shown in Figure 45. Here, the strain paths represent the corresponding deformation modes considering the material anisotropy, which are theoretically predicted from Equation 8.10 when the sheet metal is subjected to typical stress states with respect to the loading direction.   2F cos2 θS Q1 + 2Gsin2 θS Q2 − Q3 − Q4  α= (2Gcos2 θS ) Q2 + 2F sin2 θS Q1 + Q3 + Q4

(8.10a)

with Q1 = βcos2 θS + sin2 θS ,

Q2 = cos2 θS + βsin2 θS

(8.10b)

Q3 = H (1 − β) (cos 4θS + 1) , Q4 = N sin2 2θS (1 − β) where α and β denote the strain and stress paths, respectively. Figure 46(a) shows an important clue that the construction of the FFLD is essential to properly assess the global forming severity when considering the level of necking-based FLD which is about two times lower than that of the FFLD in this case. It can be also concluded that the directionality of the fracture strain at the typical stress state is not able to be adequately predicted when employing the isotropic ductile fracture criterion as shown in Figure 46(b). The strain- and stress-based fracture forming limit envelopes are represented in the space of the anisotropic stress triaxiality and the loading direction as shown in Figure 47. Since the principal stresses and strains are possible to be expressed in terms of the anisotropic stress triaxiality as discussed in chapter 4, the fracture forming limits can be evaluated in the general three-dimensional space of the principal stresses and strains with the level of anisotropic stress triaxiality. Figure 48 shows the variation of the three-dimensional fracture envelopes. In these representations, it can be seen that the level of each principal components for both cases decreases as the anisotropic stress triaxiality increases. This result is quite reasonable because an increase in the level of stress triaxiality implies an increase in the hydrostatic stress which will accelerate the growth of voids in the material.

8.2.1

Prediction of fracture strains at additional loading directions

The improvement of the model performance on the fracture prediction was confirmed in section 8.2 by comparing the fracture forming limits predicted from the proposed fracture criterion. This is quite promising because the parameters of the fracture criterion were calibrated directly from the experimental results. That is, more weight in the fracture prediction will be 57

given to the stress states and the loading directions involved in the calibration of the model parameters. It is thus necessary to validate the predictability of the fracture criterion with the experimental data not involved in the calibration of the model parameters. Additional experiments of the DP980 1.2t steel sheet have been carried out accordingly with specimens fabricated along 22.5° and 67.5° to the rolling direction. Three different types of specimen geometries were considered to induce typical stress states of the pure shear, the uniaxial tension, and the plane strain tension. Note that those specimen geometries were designed to have the same dimensions for each type considered in the previous experiments so as to disregard the possible influence of specimen geometry changes on the experimental data for the reliable comparison with the previous ones. Tensile tests were conducted three times to confirm the reproducibility from each test. With the Digital Image Correlation (DIC) method, equivalent plastic strains at the onset of fracture were measured at the material point showing the maximum equivalent plastic strain before the onset of fracture. The obtained results are summarized in Table 23. In comparison of the data evaluated from the experiments and the ones predicted from the proposed fracture criterion, it is clearly confirmed that the proposed fracture criterion is capable of predicting the onset of fracture over a wide range of stress states with great accuracy. The comparison results and the overall fracture prediction by the proposed fracture criterion are represented in Figures 49 and 50 according to the loading direction.

58

Chapter 9 Numerical Implementation 9.1

Flow rule Flow plasticity is a solid mechanics theory used to describe the plastic behavior of materials

[121]. In flow plasticity theories, the total strain in a body can be decomposed additively (or multiplicatively) into an elastic and a plastic (inelastic) parts. The elastic part of the total strain can be computed from a linear elastic or hyperelastic constitutive law: however, the determination of the plastic part of the total strain requires a flow rule and a hardening model. The plastic flow is mainly due to the rate of plastic strain increment vector and it can be assumed to be derived from a plastic potential g as: dε˜p = dλ

˜ ∂g (σ) ˜ = dλ∇g (σ) ˜ ∂σ

(9.1)

where dλ represents a positive scalar of proportionality called a plastic multiplier factor (compliance) that defines the amount of plastic straining and it can change with loading. The plastic potential determines the direction of plastic strain increment and is considered to lie in the strain space similar to the yield criterion in the stress space. Since the stress and strain variables are commonly thought of as being interchangeable, it can be allowed to make the potential function of g = 0 that describes a surface for determining the plastic flow direction in stress space. In the case that the plastic potential and the yield criterion are equal (as is normally assumed), the flow rule is termed associative. This is valid for stable materials for which, according to Drucker’s postulate [58, 59], the yield surface is convex (convexity rule) and the plastic straining in a smooth point of the yield surface takes place along the outward normal of the yield surface. This directionality of the plastic strain increment vector is referred to as a normality rule which leads to the condition of maximum plastic work [60]. The hypothesis of the maximum plastic work describes admissible directions of the plastic strain increment dεpij by assuming that the stress tensor not only must satisfy the condition of F ≤ 0, but also should be such as to maximize the increment (or rate) of plastic strain [54], i.e., ˜ : dε˜p } max {dW = σ

(9.2)

After turning a maximization problem of 9.2 into a minimization problem by changing dW to −dW , a Lagrangian function can be created by adding a Lagrangian multiplier dλ times the

59

yield criterion F : ˜ dλ) = −σ ˜ : dε˜p + dλF L (σ,

(9.3a)

˜ − h (¯ F = fy (σ) εp )

(9.3b)

˜ and dλ yields the Kuhn–Tucker conMaking L stationary with respect to variations on σ ditions as follows: ˜ ∂fy (σ) ∂L ∂F = −dε˜p + dλ = 0 → dε˜p = dλ ˜ ˜ ˜ ∂σ ∂σ ∂σ

(9.4a)

F ≤0

(9.4b)

dλ ≥ 0

(9.4c)

dλF = 0

(9.4d)

Consequently, we can then obtain the normality rule as given in Equation 9.4a. Equation 9.4d, called complementarity condition, requires that either the yield criterion is zero or dλ is zero when there is no plastic flow. The normality rule ensures uniqueness of the boundary value problem solution. Bishop and Hill [28] discussed that the normality rule is theoretically valid for polycrystalline materials. From extensive reviews of experimental yield surface results, Hecker [81] found that the normality rule is reasonable for most single phase like materials. The Associative Flow Rule (AFR) was strengthened by experimental observations of Bridgman [30, 31]. In their research, a series of tensile tests on metals was conducted in the presence of high hydrostatic pressure and the test results revealed that the pressure had no influence on the yielding of the material. In addition, Khan and Huang [97] reported that a negligible permanent volume change was shown to exist, i.e. ∆εpkk = 0. Due to the absence of pressure sensitivity in the yielding behavior, only the deviatoric stress has been involved in the formulation of the yield criterion for the materials obeying that observation. It is worth to note that the zero plastic dilatancy (zero permanent volume change) will not be violated by using the associative AFR with the employment of yield criterion without pressure sensitivity. Meanwhile, the AFR loses its validity for porous, granular, and geologic materials that exhibit pressure dependence on the yielding behavior [106, 164]: in general, those materials do exhibit a volumetric dilatation during plastic straining. For example, the volumetric dilatation predicted by the Drucker– Prager yield criterion [61] with the AFR is often somewhat larger than that evaluated from the actual experiment. It means that the inelastic volume changes in pressure-sensitive materials can be overestimated by using the AFR. Even for the metals, the AFR can be violated in the regime of substantial void growth as reported by Spitzig and Richmond [158]. In their 60

research, they showed that the AFR overestimates the plastic dilatation in the presence of superimposed hydrostatic pressure. In order to resolve the drawback of the AFR, the flow rule can be considered as non-associative as it is initially defined (non-equivalence of the yield criterion and the plastic potential).

9.1.1

Non-associative flow rule

In the framework of the AFR hypothesis, dozens of phenomenological anisotropic yield criteria have been introduced, showing the good performance in describing the material anisotropy. Recently, the inability of the AFR approach, however, has been reported by various research works [55, 142, 162, 186] in dealing with highly anisotropic materials such as AA 2008-T4, AA 2090-T3, AA 5042, and so on. In other words, the AFR concept does not always guarantee the sufficient accuracy in predicting the anisotropic material behavior even though the prediction accuracy can be also influenced by the characteristic of the constitutive model used in numerical analysis. Unlike the AFR approach, the non-associative flow rule does not enforce the artificial constraint of the equivalence between the yield criterion and the plastic potential as shown in Figure 51. In the non-AFR approach, the yield criterion defines the onset of plastic deformation while the plastic potential describes the direction of the plastic strain increment vector, respectively. From a mathematical point of view, employing two separate functions brings about an increase in the model parameters for description of the material behavior. Consequently, this approach can produce a better agreement between simulation and experimental results considering the material anisotropy. Since there is no direct link between the yield stress and the r-value in the non-AFR approach, the first and second order anomalous behaviors can be modeled even by the Hill’s 48 criterion with employing the non-AFR approach. From theses advantages and flexibility, the use of the non-AFR approach has been remarkably increasing in the field of metal plasticity [55,88,129,142,147,158,164,165,170]. In Stoughton and Yoon [164], the stability conditions for rate and temperature insensitive material models have been discussed in employing the non-AFR approach as follows: (1) To avoids a singularity for the stress rate boundary condition: H=

dh(¯ εp ) >0 d¯ εp

(the slope of the stress–strain relation)

(9.5)

(2) To avoid a singularity and ensure that the effective plastic strain rate is positive: H+

∂fy ∂fp Cpqrs >0 ∂σpq ∂σrs

61

(9.6a)

 ∂fy ∂fp Cabkl Cijcd  ∂σcd ∂σab   E˙ σ˙ ij =  C − H ep ijkl  ∂fy ∂fp  kl H+ Cpqrs ∂σpq ∂σrs     Hep ∂fy ˙ ˙ε¯p =  ∂fp ∂fy  Cijkl Ekl > 0 H+ Cpqrs ∂σij ∂σpq ∂σrs 

(9.6b)

(9.6c)

(3) To ensure that the first order plastic work rate is positive: ∂fp ˙ = σij E˙ p = σij ε¯˙p ∂fp > 0 σij > 0 → W ij ∂σij ∂σij   Hep = 0 : Pure elastic deformation p e ˙ ˙ ˙ Eij = Eij + Hep Eij →  H = 1 : Elastic–plastic deformation

(9.7a)

(9.7b)

ep

(4) To ensure that the second order palstic work rate is positive for all proportional loading:   σij = σs (t) nij σ˙ s σ˙ ij = σij → : under an arbitrary proportional loading (9.8a)  σ˙ = σ˙ (t) n σs ij

s

ij

∂fp ∂fy σ˙ kl σ˙ ij ∂σkl ∂σij ∂fy p ˙ = σij > 0 → σ˙ ij Eij = ∂σij H



σ˙ s σs

2 

∂fy σkl ∂σkl H



∂fp σij ∂σij

 >0

(9.8b)

where fy and fp denote a continuously differentiable function of the stress tensor and the plastic potential, respectively.

9.2

Explicit time integration In the previous section, the merit of the non-AFR approach has been briefly discussed

in dealing with the anisotropic material behavior. For the implementation of the non-AFR approach into finite element code, it is essential to integrate a constitutive model for the calculation of the incremental values of stress and strain considering the effect of the loading path change. In the case that the material model is not yet included in the material library of the commercial FE software, a user-defined material subroutine has to be implemented manually. The user-defined material subroutine is a sort of code written in a programming language such as FORTRAN or C. In the commercial FE software of ABAQUS, various constitutive models are already built in and it keeps updating its material library every year. ABAQUS do support both implicit and explicit user-defined material subroutines: UMAT and VUMAT. It is usually recommended to develop an explicit algorithm (VUMAT) before dealing with an implicit algorithm (UMAT) because the VUMAT does not need the calculation of the consistent tangent 62

stiffness matrix. In the UMAT, each node in the FE model depends on the neighboring nodes whereas they are independent in the VUMAT. The VUMAT includes a do loop so as to calculate the stress in every integration point in the model. ABAQUS/Explicit works on the basis of an explicit time integration scheme which uses the second order central difference scheme to integrate the equations of motion explicitly through time. In an explicit scheme, the solution at time t+∆t is only based on known quantities at time t. There are two main steps: the first step is nodal calculations and the other step is element calculations. In the first step, the accelerations at the beginning of the current increment time are evaluated based on the dynamic equilibrium equation [26]: ¨ = M−1 (Pt − It ) u Z

(9.9a)

˜ BT σdV

I=

(9.9b)



where M is the nodal mass matrix lumped for improving computational time, P is the external applied forces, I is the internal element forces, and B is the displacement differentiation matrix, respectively. After calculating the acceleration at time t, the velocities at the middle of the current increment and the displacements at the end of the increment are evaluated subsequently: ¨ t = u˙ t− ∆t + u˙ t+ ∆t = u˙ t− ∆t + ∆tt+ ∆t u 2

2

2

2

∆tt+∆t + ∆tt ¨t u 2

ut+∆t = ut + ∆tt+∆t u˙ t+ ∆t

(9.10a) (9.10b)

2

Figure 52 represents the relationship among the acceleration, velocity, and the displacement defined in the explicit time integration based on the central difference scheme. Note that the explicit time integration scheme is conditionally stable. There exist a critical size of time increment, which cannot be exceeded to obtain a stable solution. The time increment influences on results of FE analysis: a small time increment ∆t gives more stable solution while it leads to longer computational time. In ABAQUS, the critical size of time increment is automatically calculated at each explicit time integration by the length of the smallest element size Le and the velocity of the dilatational wave speed cd which depends on a Young’s modulus E, Poisson’s ratio ν, and specific mass density ρ, e.g. cd = [E/(ρ − ρν 2 )]−1/2 for the shell element. The critical size of time increment ∆t is then defined as:  ∆tcr = min

63

Le cd

 (9.11)

It is worth to highlight that the time step increment in the explicit time integration must be smaller than the critical value for the stable solution, i.e. ∆t ≤ ∆tcr . In the second step, ABAQUS subsequently provides a stress state, a strain increment and other information corresponding to the previously converged step from FE software for each element to the userdefined material subroutine. Here, the strain increment, given in Equation 9.12b, is computed by using the time increment evaluated during the process of the explicit time integration: T

ε˙t+ ∆t = B u˙ t+ ∆t 2

2

  ∆tt+∆t + ∆tt ¨t u˙ t− ∆t + =B u 2 2 T

(9.12a)

∆εt+ ∆t = ε˙t+ ∆t ∆tt+∆t 2

(9.12b)

2

In the user-defined material subroutine, the stress state is then updated from the numerical integration of the constitutive equation based on an imposed strain increment:   σt+∆t = σt + ∆σ = σt + ∆tt+∆t σ˙ t+ ∆t ← σ˙ t+ ∆t = F ε˙t+ ∆t 2

2

(9.13)

2

The assemble nodal internal forces at time step t + ∆t can be then computed based on the updated stress state: Z

˜ t+∆t dV BT σ

It+∆t =

(9.14)



Consequently, these two main steps are repeated to update the strain and stress states for each element so as to describe the material behavior during the explicit FE analysis.

9.2.1

Stress update algorithm

In the previous section, two numerical steps are reviewed about the process for the explicit FE analysis in ABAQUS. This section will discuss the stress update algorithm which features the non-AFR approach. There exist various numerical techniques for the determination of the plastic multiplier factor ∆λ during the stress update process: A backward Euler method (a fully implicit approach); A forward Euler method (a fully explicit approach); A tangent cutting plane method (a semi-explicit approach) [156]; Next Increment Corrects Error (NICE) method (a fully explicit approach) [79] and etc. Among those methods, the backward Euler method is the most favorably used one in the process of the stress integration due to its unconditional stability and quadratic convergence rate which is mainly related to the inherent characteristic of the Newton–Raphson iteration method. This method, however, requires second order derivatives of the yield criterion and the plastic potential, which gives rise to difficulties in the formulation of the stress update algorithm, particularly in the case that the constitutive equations used 64

become complicated. It is worth to note that explicit integration approaches can give acceptable results when those are used with very small time increments in the explicit time integration. In other words, the degree of violation of the consistency condition will be elevated if large time increments are used. Cardoso and Yoon [39] reported that the forward Euler method is possible to fulfill the same degree of accuracy as can be achieved by the backward Euler method when used in the explicit time integration with small time increments. Using very small time increments in the explicit time integration, however, can generate some disadvantages in computational efficiency in the case that the FE analysis becomes complicated. As a compromise between the computational time and the numerical accuracy, the tangent cutting plane method can be considered as an alternative integration method which is a sort of the semi-implicit approaches. In this method, there is no need to obtain the second order derivatives of the yield criterion and the plastic potential because of its explicit characteristic in regard to the plastic multiplier factor ∆λ. The simplicity of this stress integration method, in general, limits the user towards using considerably smaller time increments so as to reduce the violation of the consistency condition. The tangent cutting plane method enforces the normality at the beginning of each increment. The consistency condition is satisfied at the end of the step so as to avoid the drift of stress from the yield surface. The consistency condition can be derived by the linearization of the yield criterion around the actual stress state. Recall the definition of the yield criterion F : ˜ ε¯p ) = fy (σ) ˜ − h (¯ F (σ, εp )

(9.15)

The consistency condition of the equation above can be then expressed as: ˜ + ∆σ, ˜ ε¯p + ∆¯ ˜ + ∆σ) ˜ − h (¯ ∆F = F (σ εp ) = fy (σ εp + ∆¯ εp ) = 0

(9.16)

Based on a total strain increment and the non-AFR approach, the stress and plastic strain increments can be defined as: ˜ = C : (∆ε˜ − ∆ε˜p ) ∆σ ∆ε˜p = ∆λ

∂fp ˜ ∂σ

(9.17a) (9.17b)

where fp is the plastic potential. In the non-AFR approach, the relationship between the plastic multiplier factor ∆λ and the equivalent plastic strain increment ∆¯ εp can be derived based on the principle of plastic work equivalence: ˜ ∆¯ ˜ : ∆˜ fy (σ) εp = σ εp → ∆¯ εp = 65

˜ ∂fp (σ) ˜ fp (σ) ˜ ∂σ = ∆λ ˜ ˜ fy (σ) fy (σ)

˜ : ∆λ σ

(9.18)

˜ ∂fp (σ) ˜ = fp (σ) ˜ ∂σ ˜ T (= C∆ε) ˜ yields the based on the Euler’s theorem. Linearization at the trial stress state σ ˜ : where fp is a first order homogenous function that obeys the relation of σ

˜ in Equation 9.17a with respect to the increment of the equivalent plastic stress difference dσ strain: i.e., ˜ = −C : ∆λ ∆σ

∂fp ∂σ ˜

(9.19)

The consistency condition linearized at the trial stress can be expressed as: ˜ + ∆σ) ˜ − h (¯ ˜ − h (¯ ∆F = fy (σ εp + ∆¯ εp ) = fy (σ) εp ) +

∂fy dh ˜ − p ∆¯ εp : ∆σ ˜ ∂σ d¯ ε

∂fy ∂fp dh ˜ − h (¯ = fy (σ) ε )− : C : ∆λ − p ∆¯ εp = 0 ˜ ˜ ∂σ ∂σ d¯ ε

(9.20)

p

Substituting Equation 9.18 into Equation 9.20 eventually gives the increment of the plastic multiplier factor dλ as: ∆λ =

˜ − h (¯ ˜ − h (¯ fy (σ) εp ) fy (σ) εp ) = ∂fy fp ∂fp dh fp ˜ :C:n ˜ +H m :C: + p ˜ ˜ ∂σ ∂σ d¯ ε fy fy

(9.21)

Using the equation above, the stress state is eventually updated for a given strain increment ˜ The numerical algorithm of the tangent cutting plane method is summarized in Table 24 ∆ε. and its graphical representation is given in Figure 53.

66

Chapter 10 Experimental Validation This chapter discusses the fracture predictability of the proposed anisotropic ductile fracture criterion for a structural level. The validity of the model performance will be verified by comparison of the experimental and numerical simulation results of the square cup drawing test. Each fracture forming limit proposed is employed in order to not only predict the onset of fracture and the crack location observed in the experiment but also confirm the model performance under the non-proportional stress state.

10.1

Square cup drawing test

In order to validate the fracture predictability of the proposed fracture criterion, the square cup drawing test was conducted with the test condition given in Table 25. All the drawing tests were designed to make the specimen being fractured during the drawing process. In the tests, 3,000 kN hydraulic servo press was used. Specifications of the test machine are summarized in Table 26. Figure 54 shows the tested specimen with a crack propagated at the punch edge of the sheet. The non-symmetric crack propagation was observed, which may result in non-uniform distribution of blank holding force. The punch stroke to fracture was evaluated as around 21 mm from the tests.

10.2

Square cup drawing simulation

All the FE simulations were performed in the environment of ABAQUS/Explicit V6.10. The crack initiation and propagation are simulated using the element deletion technique. The element is deleted when the damage index predicted by the proposed fracture criterion reaches unity at each integration point. The punch, die, and blank-holder were modeled as the discrete rigid body in ABAQUS, with the mesh size of 2mm x 2mm. In order to describe the nonsymmetric crack propagation, simulations were performed using two different FE models of the square sheet: a quarter and a half FE models. Schematic descriptions of each square cup drawing simulation are represented in Figure 55. In the case of a half model, a bias of 5 mm was imposed as illustrated in Figure 55(b) to generate non-uniform blank-holding force during the drawing simulation. The geometries of both FE models were discretized using three-dimensional brick elements with reduced integration (C3D8R in ABAQUS element library). Each FE model 67

was made with various mesh size so as to confirm the dependence of mesh size on the fracture initiation. Six element layers were made through thickness direction for all cases to compensate the computational efficiency. The total number of elements and nodes are summarized in Table 27 for each case. The friction coefficient of 0.1 was assigned in the FE simulations and assumed the same everywhere between the contact regions of the FE models. In the FE simulations, the non-associated flow rule was adopted to incorporate the effect of a different directionality of the r-value and the yield stresses without losing advantages of quadratic functions. In this approach, both plastic potential and yield criterion were formulated to have the same form of the Hill’s 48 criterion. For the stress update process, the constructed loci are represented in Figure 56. The tangent cutting plane method was employed in consideration of the computational efficiency as well as the numerical accuracy. The material behavior was predicted by the hardening curves evaluated from the combination of the tensile and hydraulic bulge tests. For confirmation of whether the stress–strain response in the large deformation range is reliable, validation was performed by comparing the load–displacement curves evaluated from the actual test and simulation: three different tests including pure shear, uniaxial tension, and the plane strain tension were used for validation. Two simulation groups termed as (A) and (B) were introduced according to the flow stress curve used in FE analysis as represented in Figure 57. Following the recommendation of Lou and Huh [118] and Dunand and Mohr [63], a very fine mesh of eight-node 3D linear brick elements with reduced integration was built up around the critical area of necking for each specimen. An element edge length of about 0.1 mm was assigned and eight element layers were made through the half-thickness. For the computational efficiency, symmetric conditions were considered in the FE modeling as shown in Figure 58. From the simulation results, as can be seen in Figure 59, it was clearly confirmed that the hardening behavior predicted from the combination of the tensile and hydraulic bulge tests is able to describe the material deformation behavior with great accuracy whereas the flow stress curve evaluated from tensile tests only overestimates the loads evaluated for each case. The results of FE analysis also showed that the fracture strokes predicted from the FE simulation were lower than those evaluated from the actual fracture tests. The reason behind this is that the fracture strains used for the calibration of the model parameters are obtained from the surface of a specimen such that the onset of fracture is somewhat overestimated in the FE simulation.

68

10.3

Comparison of the experimental and simulation results

All square cup drawing simulations were performed using the flow stress curve obtained from the combination of the tensile and hydraulic bulge tests as its validity was verified in subsection (p)

10.2. A projection method proposed by Stoughton and Yoon [166], i.e. σij = σij h (¯ εp ) /fy (σij ), was taken into consideration to compensate possible numerical violation of the stress–strain relationship defined by the constitutive model during the space–time integration approximation. Figures 60 and 61 show the simulation results for both simulation groups: one using a quarter FE model; and the other using a half FE model. As expected, particularly for a half FE model with a bias of 5 mm, the non-symmetric crack propagation was predicted although a quarter FE model showed symmetric crack propagation as represented in Figure 62. Simulation results of both cases demonstrate that the punch stroke to fracture decreases as the mesh size decreases. This tendency is quite reasonable because the decrease in the mesh size of the FE model implies the increase in the number of numerical integration points. In other words, numerical analysis using the FE model with the coarse mesh demonstrates the deformation behavior in an average sense. Figure 63 shows distributions of the localized equivalent plastic strain through the thickness with the variation of the mesh size along the punch radius. Meanwhile, both simulation groups predicted the fracture location in the blank considering the material orientation appropriately. This result can reveal the fact that the proposed fracture criterion has a sufficient performance in predicting the fracture initiation for a certain stress state that the material mainly undergoes in consideration of the material anisotropy. In order to evaluate the dependence of the mesh size on the fracture initiation, the relationship between the punch stroke to fracture and the mesh size normalized by the punch radius was investigated because the fracture takes place at the location near the edge of punch radius during the drawing process. The tendency of the punch stroke to fracture was described for each FE group using an exponential function given in Equation 10.1 as shown in Figure 64.  Pf = P0 + P1 exp

x P2

 (10.1)

where x denotes the global mesh size normalized by the punch radius. P0 , P1 , and P2 are the model parameters which are given in Table 28 for both cases. From Equation 10.1, the punch stroke to fracture can be evaluated when the mesh size becomes zero. If we assume that the total number of elements is infinite, the range of punch stroke to fracture is from 18.92 mm (for a quarter model) to 21.80 mm (for a half model). For the case of a half model with the mesh size of 0.5 mm, the fracture evaluation was performed by using the strain- and 69

stress-based FFLD as well as the PEPS FFLD. Each evaluation is represented in Figures 65–71 and the results demonstrate that the fracture initiation is necessary to be evaluated using the stress-based FFLD and the PEPS FFLD particulary for the non-porportional loading state. Comparison of the experimental and numerical simulation results gives a clear conclusion that the proposed fracture criterion has a sufficient capability to predict the fracture initiation over a wide range of stress states in consideration of the material anisotropy and non-directionality of the equi-biaxial fracture strain.

70

Chapter 11 Conclusions and Future Works 11.1

Conclusions

Anisotropy effect on the fracture strain was investigated for the DP980 1.2t steel sheet. Tensile tests were conducted at the three loading directions of (0°, 45°, 90°) in consideration of typical loading conditions including the pure shear, the uniaxial tension, and the plane strain tension. From the tests results, it was clearly found that the fracture strain varies according to the loading direction. The non-consistent tendency was observed between the fracture strains for typical loading state according to the loading direction. With use of tests results, the quantitative evaluation on the fracture strain was carried out in accordance with the maximum principal stress direction not only by developing the anisotropic ductile fracture criterion but also by introducing forming limit criteria for the anisotropic material as follows: (1) A new uncoupled ductile fracture criterion is developed to provide a phenomenological model for considering non-directionality of the equi-biaxial fracture strain as well as the effect of anisotropy at the macroscopic level particularly for sheet metal forming application. In developing the anisotropic ductile fracture criterion, the Hill’s 48 criterion was employed to take account of the loading history and the anisotropic effect on the equivalent plastic strain at the onset of fracture accurately. The principal stresses were expressed in terms of (ηH , LP , σ ¯H ) to include the effect of directionality for typical loading states on the material orientation. From the combination of the stress triaxiality defined by the Hill’s 48 criterion and the Lode parameter, it was shown that the corresponding principal stresses are uniquely determined in Haigh–Westergaard space for the anisotropic material. (2) In addition to the construction of the anisotropic ductile fracture criterion, the concept of the generic anisotropic forming limit on the strain-based FFLD was employed to propose a representative fracture-based forming limit for anisotropic materials. The constructed strain-based fracture envelope was shown to be capable of describing ductile behavior under the proportional loading condition. (3) For a direct application of the fracture prediction to sheet metal forming processes, the strain-based FFLD was theoretically derived from the proposed fracture criterion by assuming the proportional loading and plane stress conditions. An attempt was made to consider the complicated path dependence in evaluating the forming severity of the ma71

terial, the stress-based fracture envelope was introduced in the space of the major and minor limit stresses with the same assumption applied in the derivation of the strain-based FFLD. (4) For easy application of the forming limit criterion to the automotive industry, the strainbased FFLD was converted into the PEPS FFLD which is the transformation of the path dependent strain-based FFLD to the path independent space similar to the path independence of the PEPS FLD. (5) All fracture-based forming limit criteria transformed from the anisotropic ductile fracture criterion proposed were capable of describing the cut-off value for the stress triaxiality as well as predicting the fracture initiation, especially for the shear ductile fracture of metals over a wide range of loading conditions from compressive loading, and shear to equi-biaxial tension. (6) Experimental validations were performed in the viewpoint of a structural and a specimen level through the square cup drawing and the tensile test at additional loading direction. In comparison of the experimental results and the ones predicted from the proposed fracture criterion, it was clearly shown that the proposed fracture criterion has a sufficient capability in describing the instance of fracture initiation over a wide range of stress states in consideration of the material anisotropy and non-directionality of the equi-biaxial fracture strain. (7) Comparison of the strain-based FFLD and the PEPS FFLD revealed that the PEPS FFLD can predict the fracture initiation regardless of the strain path change whereas the performance of the strain-based FFLD can be valid only for the case that the material point is under the proportional loading condition. Consequently, the PEPS FFLD can play a key role in predicting the fracture initiation during complicated forming processes of AHSS sheets.

11.2

Future works

(1) A microscopic investigation should be carried out to explain the reason for the dependence of the fracture strain on the sheet metal orientation. (2) The effect of strain rate and temperature on the fracture strain should be investigated considering the material anisotropy and incorporated into the proposed fracture criterion.

72

Appendix A. Generalization of the Anisotropic Stress Triaxiality The anisotropic stress triaxiality suggested in section 4.3 can be generalized for application to the aluminum alloy based on the non-quadratic yield criterion. In the generalization of the anisotropic stress triaxiality, it is worth to mention that the non-quadratic anisotropic yield criterion employed have to deal with the full stress states that the material usually undergoes. This requirement is essential to transform the principal stress states into the space of the anisotropic stress triaxiality, the Lode parameter, and the anisotropic yield criterion employed.

A.1 Anisotropic stress triaxiality defined by the Yld2004-18p criterion Barlat et al. [19] proposed a generalization of the Hosford criterion for a pressure-independent material under general stress states with two linear transformations:

φ

0 σαβ



1,3  X ˜0 ˜00 m 0 ˜00 m ˜ = φ Si , Sj = σA Si − Sj = 4¯



(A.1)

i,j

where m is a constant coefficient mainly associated with the crystal structure and the subscripts α and β stand for x, y, and z. Note that the curvature of the rounded vertices of the yield surface increases, as the m increases [187]. Here, a reference frame (x, y, z) is assumed to be attached to the material symmetry axes as shown in Figure 17. In Equation A.1, S˜i0 and S˜j00 are the principal values of the two transformed stress deviator ˜s0 and ˜s00 which can be written in a matrix form as   s˜xx       s˜yy       s˜zz ˜s = Cσ ˜0 =   s˜yz       s˜xz      s˜xy

             



0

−c12 −c13

0

0

0

  0 −c23 0 0 0  −c21    −c31 −c32 0 0 0 0 =    0  0 0 c44 0 0         0  0 0 0 c55 0       0 0 0 0 0 c66

73

                           

 0  σxx      0 σyy       0 σzz  0   σyz     0  σxz       0 σxy

(A.2)

with the relevant symbols (prime and double prime) for each transformation: i.e. c0ij for ˜s0 and c00ij for ˜s00 , respectively. The transformation tensor can also apply on the principal stress tensor of the Cauchy stress state as ˜s = CTσ ˜ = Lσ ˜ ⇒ ˜s = Lσ ˜ = LT∗ σ ˜ p = L∗ σ ˜p

(A.3a)

where       1 T=  3     



2 tx2 + c

−1 −1 0 0 0

2



  −1 0 0 0    −1 −1 2 0 0 0    0 0 0 3 0 0    0 0 0 0 3 0   0 0 0 0 0 3

−1

2

(txy + sz)2 2 ty 2 + c

(txz − sy)2

  (txy − sz)2 (tyz + sx)2     2 2 2+c 2  (txz + sy) (tyz − sx) tz T∗ =     2 + c (tyz − sx) 2+c  (txy − sz) (txz + sy) ty (tyz + sx) tz      tx2 + c (txz + sy) (txy + sz) (tyz − sx) (txz − sy) tz 2 + c    tx2 + c (txy − sz) (txy + sz) ty 2 + c (txz − sy) (tyz + sx)

(A.3b)

             

(A.3c)

ˆ is a unit ˜ p = { σ1 σ2 σ3 }T , t = 1 − cos θr , c = cos θr , s = sin θr , u where σ ˆ = xˆi + yˆ j + zk vector for a rotation axis and θr is the rotation angle as shown in Figure 17. From Equations A.2 and A.3a, the following relation can be obtained:   s˜xx       s˜yy       s˜zz   s˜yz       s˜xz      s˜xy

             

  c12 (σxx − σyy ) + c12 (σzz − σyy ) + c13 (σxx − σzz ) + c13 (σyy − σzz )       c21 (σyy − σxx ) + c21 (σzz − σxx ) + c23 (σxx − σzz ) + c23 (σyy − σzz )      1  c31 (σyy − σxx ) + c31 (σzz − σxx ) + c32 (σxx − σyy ) + c32 (σzz − σyy ) =  3   3c44 σyz             3c55 σxz           3c66 σxy

             

(A.4)

            

The Cauchy stresses in Equation A.4 can be expressed in terms of the principal stresses as discussed in section 4.3. Equation A.4 can be thus rewritten as

74

  s˜xx       s˜yy       s˜zz   s˜yz       s˜xz      s˜xy

  c12 h + (c13 − c12 ) f + c13 (−g)       c21 (−h) + (c21 − c23 ) g + c23 f      σ1 − σ3  (c32 − c31 ) h + c31 g + c32 (−f ) =   3   3c44 l             3c55 m           3c66 n              

             

(A.5)

            

where the subfunctions of f, g, h, l, m, and n are rigorously derived by the tensor transformation in aid of the Rodrigues’ rotation formula as shown in Equation 4.32c. The principal values of ˜s, ˜ are the roots of the characteristics equation: which defines the diagonal tensor representation S, P (S˜k ) = S˜k3 − H1 S˜k2 − H2 S˜k − H3 = 0

(A.6a)

where H1 , H2 , and H3 are the first, second and third invariants defined as H1 = s˜kk = H2 =

s˜xx + s˜yy + s˜zz 3

(A.6b)

1 (˜ sii s˜jj − s˜ij s˜ji ) = s˜xx s˜yy + s˜yy s˜zz + s˜zz s˜xx − s˜2xy − s˜2yz − s˜2zx 2

H3 = det(˜ sij ) = s˜xx s˜yy s˜zz + 2˜ sxy s˜yz s˜zx − s˜xx s˜2yz − s˜yy s˜2zx − s˜zz s˜2xy

(A.6c) (A.6d)

Solving the above characteristic equation gives the principal values of ˜s in terms of its invariants:       ˜   1 S    1    H  1 ˜ = S2 1   3         1  S˜  3

    

+2

p

−Q

              

   

 θA cos 3   θA − 2π cos 3   θA + 2π cos 3 

       

(A.7a)

      

where 2H13 − 9H1 H2 + 27H3 3H2 − H12 , R= , θA = cos−1 Q= 9 54

R

!

p −Q3

(A.7b)

Because 0 ≤ θA ≤ π, the principal values are ordered, S˜1 ≥ S˜2 ≥ S˜3 . Substituting Equation A.5 into Equation A.6 yields   σ1 − σ3  c12 h + (c13 − c12 ) f + c13 (−g) + c21 (−h) + (c21 − c23 ) g + c23 f  H1 = 3 +c g + (c − c ) (−h) + c (−f ) 31

31

75

32

32

(A.8a)

    c12 h + (c13 − c12 ) f + c13 (−g) c21 (−h) + (c21 − c23 ) g + c13 f     3 3         c (−h) + (c − c ) g + c f c g + (c − c  21 21 23 13 31 31 32 ) (−h) + c32 (−f )   + 3 3 H2 = (σ1 − σ3 )2      c h + (c − c12 ) f + c13 (−g) c g + (c − c ) (−h) + c (−f ) 12 13 31 31 32 32   +   3 3       −c266 n2 − c244 l2 − c255 m2

            

(A.8b)

           

     c31 g + (c31 − c32 ) (−h) + c32 (−f ) c12 h + (c13 − c12 ) f + c13 (−g) c21 (−h) + (c21 − c23 ) g + c13 f     3 3 3       c12 h + (c13 − c12 ) f + c13 (−g) 2 2 3 2 2 2 2 2 2 +2c66 n c44 l c55 m − H3 = (σ1 − σ3 ) c44 l  3          c31 g + (c31 − c32 ) (−h) + c32 (−f ) 2 2 c21 (−h) + (c21 − c23 ) g + c13 f 2 2   c55 m − c66 n − 3 3

        

(A.8c)

       

By plugging Equation A.8 into Equation A.7b, the sub-quantities of Q, R, and θA is simply expressed as 2H13 − 9H1 H2 + 27H3 3H2 − H12 = (σ1 − σ3 )2 A, R = = (σ1 − σ3 )3 B 9 54       (σ1 − σ3 )3 B R B −1 −1 −1 q  = cos p √ θA = cos = cos  −Q3 −A3 −(σ1 − σ3 )6 A3

Q=

(A.9)

where A and B are the non-dimensional quantities associated with the normalized components of the rotation axis u ˆ = (x, y, z), the Lode parameter LP , and the rotation angle θr . With the transformation procedure introduced, the Yld2004-18p criterion, defined by the principal values of the transformed stress deviators of ˜s0 and ˜s00 , is expressed in terms of the principal stresses of the Cauchy stress tensor accordingly:  0   00  1,3 X √ √ θA + 2π (i − 2) θA + 2π (j − 2) B 0 − B 00 m 0 00 φ = |σ1 − σ3 | − 2 −A cos + 2 −A cos 3 3 3 (A.10) i,j m

m = 4¯ σA

with the symbols (prime and double prime) for each transformation of c0ij and c00ij . It is worth to mention that the equivalent stress of the Yld2004-18p criterion is unequivocally associated with the diameter of the Mohr’s circle. This is sound from the point of view that the Yld2004-18p criterion starts from a framework of the maximum shear stress theory. The anisotropic stress triaxiality based on the Yld2004-18p criterion is defined by the principal stresses of the Cauchy stress tensor: ηA =

=

σm = σ ¯A

σ1 + σ2 + σ3 1 )   1 ( 1,3  0   00  0 − B 00 m m √ θ + 2π (i − 2) θ + 2π (j − 2) 1 m P √ B A A 0 3 |σ1 − σ3 | − 2 −A00 cos + 2 −A cos 4 3 3 3 i,j

σ1 + σ2 + σ3 3 |σ1 − σ3 | Y

76

(A.11)

A unique relation between (σ1 , σ2 , σ3 ) and (ηA , LP , σ ¯A ) can be easily derived similary to analogy in Equations 4.34–4.36. Based on the definition of the Lode parameter and the equivalent stress of the Yld2004-18p, following relations can be obtained:

σ1 + σ3 = 2σm −

LP (σ1 − σ3 ) 2σm LP (σ1 − σ3 ) LP σ1 + σ3 = − = 2ηA − → 3 σ ¯A σ ¯A 3¯ σA 3Y

σ ¯A = |σ1 − σ3 | Y →

1 |σ1 − σ3 | = σ ¯A Y

(A.12a)

(A.12b)

where 1 m )   1 ( 1,3  0   00  0 00 √ m θA + 2π (i − 2) θ + 2π (j − 2) 1 m P √ B −B 0 (A.12c) Y = − 2 −A00 cos A + 2 −A cos 4 3 3 3 i,j

After simple mathematical manipulations, the principal stresses can be defined as  3 − LP ηA + σ ¯A 6Y   2LP 2LP σ ¯A = ηA + σ ¯A σ2 = σm + s2 = σm + 6Y 6Y   3 + LP 3 + LP σ3 = σm + s3 = σm − σ ¯A = ηA − σ ¯A 6Y 6Y

3 − LP σ1 = σm + s1 = σm + σ ¯A = 6Y



(A.13a)

(A.13b) (A.13c)

where 1 )   1 ( 1,3  0   00  0 − B 00 m m √ θ + 2π (i − 2) θ + 2π (j − 2) 1 m P √ B A 0 Y = − 2 −A00 cos A + (A.13d) 2 −A cos 4 3 3 3 i,j

From Equation A.13, yield loci lying on three symmetric planes can be easily described concerning the maximum principal stress direction. With an assumption of the plane stress condition, the normalized yield stress ratios are expressed in terms of the Lode parameter: σ1 3 − LP 1 = ηA + = , σ ¯A 6Y Y σ1 3 − LP 1 − LP σ2 = 0 → = ηA + = , σ ¯A 6Y 2Y σ2 2LP 1 − LP σ1 = 0 → = ηA + =− , σ ¯A 6Y 2Y σ3 = 0 →

σ2 2LP 1 + LP = ηA + = σ ¯A 6Y 2Y σ3 3 + LP 1 + LP = ηA − =− σ ¯A 6Y 2Y σ3 3 + LP 1 = ηA − =− σ ¯A 6Y Y

(A.14a) (A.14b) (A.14c)

where Y is a non-dimensional quantity associated with the normalized components of the rotation axis u ˆ = (x, y, z), the Lode parameter LP , and the rotation angle θr . Note that the 77

anisotropic stress triaxiality generalized by the Yld2004-18p criterion can describe the various anisotropic behavior of materials including an isotropic behavior by imposing certain conditions to the material constants, which are briefly summarized in Figure 73. The validity of the derived formula was confirmed by constructing yield loci of AA6111-T4 whose material parameters for the Yld2004-18p criterion were obtained by Barlat et al. [19]. Figure 74 shows the yield loci constructed under the plane stress condition.

A.2 Construction of the uncoupled anisotropic fracture criterion based on the Yld91 criterion In an attempt to widen the applicability of the anisotropic stress triaxiality to various materials including aluminum alloys, an uncoupled anisotropic fracture criterion is developed by employing the anisotropic indicator defined by the Yld91 criterion since the anisotropic indicator generalized by the Yld2004-18p criterion is quite complex to be readily applied in general cases. The anisotropic indicator defined by the Yld91 criterion is expressed as σm 1 σ1 + σ2 + σ3 = r m σ ¯A 3 m |S2 − S3 | + |S3 − S1 |m + |S1 − S2 |m 2

ηA =

(A.15a)

where m 2¯ σA = |S2 − S3 |m + |S3 − S1 |m + |S1 − S2 |m

(A.15b)

Here, S1 , S2 , and S3 are the principal values of the Isotropic Plastic Equivalent (IPE) transformed stress tensor s, and the exponent m denotes an isotropic parameter that can assume any positive and real value greater than unity. In the Voigt notation, the IPE transformed stress vector ˜ s is defined   sxx       syy       szz ˜s = Lσ ˜=   syz       szx      sxy

as              



      =                   

(c2 + c3 )/3

−c3 /3

−c2 /3

0

0

0

−c3 /3

(c3 + c1 )/3

−c1 /3

0

0

0

−c2 /3

−c1 /3

(c1 + c2 )/3

0

0

0

0

0

0

c4

0

0

0

0

0

0

c5

0

0

0

0

0

0

c6

                           

 σxx       σyy       σzz  (A.16)  σyz       σzx       σxy

where L is the linear transformation operator proposed by Barlat et al. [21] and ci (i = 1, . . . , 6) are the parameters describing the anisotropy of the metallic material. In aid of tensor transformation based on the Rodrigues’ rotation formula [151] with the Lode parameter, the IPE transformed stress tensor s can be written as 78

   s = (σ1 − σ3 )   

c3 h − c2 g 3 c6 n c5 m

c6 n c1 f − c3 h 3 c4 l

c5 m





sxx sxy sxz

    = s c4 l   xy c2 g − c1 f  sxz 3

syy syz



  syz   szz

(A.17)

where f , g, h, l, m, and n are the sub-functions in Equation 4.32c. The principal values of the IPE transformed tensor s can be expressed by its invariants as    θA     cos       S  3      1    p θA − 2π cos = 2 −Q S2    3          S   θ + 2π  A 3   cos 3

       

(A.18a)

      

where

θA = cos−1

R

! , Q=

p −Q3

H3 H2 , R= 3 2

     −c3 h + c1 f −c3 h + c1 f −c1 f + c2 g −c2 g + c3 h + 3 3 3 3    −c2 g + c3 h −c1 f + c2 g + − c26 n2 − c24 l2 − c25 m2 3 3



   −c3 h + c1 f −c1 f + c2 g −c2 g + c3 h + 2c26 n2 c24 l2 c25 m2 3 3 3       −c3 h + c1 f −c1 f + c2 g 2 2 −c2 g + c3 h 2 2 2 2 c4 l − c5 m − c6 n − 3 3 3



    H2 = (σ1 − σ3 )2      H3 = (σ1 − σ3 )3  

(A.18b)



   

   

(A.18c)

(A.18d)

where H2 and H3 are the second, and the third invariants of the IPE transformed stress tensor s, respectively. Note that the first invariant H1 of the IPE transformed stress tensor s is zero (skk = 0) and the principal values are in the order of S1 ≥ S2 ≥ S3 because 0 ≤ θA ≤ π. Substituting Equation A.18 into Equation A.15 gives σ1 + σ2 + σ3 1   −1      m     1 √ 1 m π − θA m θA θA + π m m |σ1 − σ3 | −3Q sin + sin + sin 2 3 3 3 1 σ1 + σ2 + σ3 = 3 |σ1 − σ3 | Y

ηA =

σm 1 = σ ¯A 3

(A.19)

It is worth to mention that the Yld91 criterion can reduce to special cases such as the Tresca [173], the von Mises [177], the Hosford [90], and the Hill’s 48 [83] criteria. For example, when all components in the linear transformation operator L are set to unity, the IPE transformed tensor s reduces to the stress deviator tensor σ 0 . In this case, the anisotropic stress triaxiality 79

proposed can reduce to the isotropic stress triaxiality by setting the exponent of m to 2. In addition, the following explicit relation between the anisotropic parameters of the Yld91 and the Hill’s 48 criteria can be obtained when the value of the exponent m is set to 2 [146]: 2c21 + c1 c2 + c1 c3 − c2 c3 3 , L = c24 6 2 2c22 + c1 c2 + c2 c3 − c1 c3 3 G= , M = c25 6 2 2 3 2c + c1 c3 + c2 c3 − c1 c2 , N = c26 H= 3 6 2 F =

(A.20)

The generalized anisotropic stress triaxiality based on the Yld91 criterion, therefore, has notable applicability for both isotropic and anisotropic cases. For the characterization of the stress states, it is now straightforward to make an explicit relation between (σ1 , σ2 , σ3 ) and (ηA , LP , σ ¯A ) because the proposed anisotropic stress triaxiality and the Lode parameter are expressed by the principal stresses of the Cauchy stress tensor σ together with the fact that the 0 = 0): first invariant of the stress deviator tensor σ 0 is zero (σkk

        σm + σ10 σ    1     = σ2 σm + σ20         σ     σ + σ0 3 m 3

 3 − LP   ηA +   6Y    2LP = ηA +   6Y        3 + LP   η − A 6Y     

       

σ ¯A

(A.21)

      

Here, the principal stresses are with the order of σ1 ≥ σ2 ≥ σ3 . From the characterization of the stress states with the anisotropic stress triaxiality, it is possible to describe yield loci lying on three symmetric planes whose axes are coincident with the material symmetry axes as shown in Figure 17. For example, when the plane stress condition is assumed as σ3 = 0 concerning sheet metal forming application, the unit vector of the rotation axis is expressed as ˆ = k, ˆ such that the sub-functions in Equation 4.32c reduce to u ˆ = xˆi + yˆj + z k      LP − 1 1 + LP LP − 1 2 2 cos θS , g = cos θS − , f =1+ 2 2 2     1 − LP LP − 1 h= cos 2θS , n = sin 2θS , l = m = 0 2 4 

(A.22)

In this case, the rotation angle θr strands for the loading direction θS as shown in Figure 75. The principal stresses lying on the normal surface of the sheet metal are then defined as    σ1  =  σ  2

 

  1 1 + L P    1 −1     0 m   0   1  0 m √ 2 π − θA θA θA + π m m 1 m 0 −3Q sin + sin + sin 2 3 3 3 σ ¯A

80

(A.23a)

where   1 −1   m   m   m  1 √ 0 0 0 1 π − θ θ θ + π A A A m m sin Y 0 = −3Q0 + sin + sin (A.23b) 2 3 3 3       LP − 1 (P1 + P2 )2 P2 LP − 1 2 1 2 2 cos 2θS − − P3 P1 + c3 − c6 sin 2θS (A.23c) 3 2 27 12 2 " ! #  2 0 1 L − 1 P + P 1 R 1 2 P R0 = c2 sin2 2θS − (P2 − 3P3 ) (P1 + 3P3 ) , θA0 = cos−1 p (A.23d) 3 8 6 2 18 −Q0 3 Q0

1 = 9



 P1 = c1

LP − 1 2



      1 − LP LP − 1 1 + LP cos2 θS + 1 , P2 = c2 cos2 θS + , P3 = c3 cos 2θS (A.23e) 2 2 6

Concerning the equi-biaxial stress state of σ1 = σ2 = σb together with σ ¯A = σb , Y 0 in Equation A.23 becomes unity regardless of the values of the anisotropic parameters ci (i = 1, . . . , 6) and the exponent m as well as the loading direction θS . This makes the anisotropic stress triaxiality become a constant value of 2/3. The loss of the directionality on the equibiaxial stress state can be reviewed from the definition of the anisotropic stress triaxiality as well. As an anisotropic indicator, the anisotropic stress triaxiality can be expressed as a quantity scaled from the isotropic stress triaxiality as follows: ηA =

σm σm σ σ ¯v ¯v = = ηv σ ¯A σ ¯v σ ¯A σ ¯A

(A.24)

where σ ¯v and σ ¯A represent the equivalent stresses of the von Mises isotropic yield criterion and the anisotropic yield criterion, respectively. A scaling factor defined by the yield stress ratio of σ ¯v /¯ σA becomes unity when the equivalent stress of σ ¯A is set to the equi-biaxial yield stress σb , which leads to ηA = ηv = 2/3. Except for the equi-biaxial stress state, the scaling values depend on the Lode parameter LP , the anisotropic parameters ci (i = 1, . . . , 6), and the exponent m:

ηA = ηv

q L2P + 3   1 −2     0 m   0 m  1 0 m √ π − θ θ θ + π 1 m m A −3Q0 sin + sin A + sin A 2 3 3 3

(A.25)

Figure 76 shows the distribution of the scaling values in general stress states according to three loading directions of 0°, 45°, and 90° from the rolling direction of the sheet metal. Here, the parameters of the Yld91 criterion are calibrated using the ratios of the normalized yield stresses given in Table 21 and are summarized in Table 29. It can be simply confirmed that the directionality of the equi-biaxial stress state (LP = 1) vanishes for the reason that the scaling value becomes unity regardless of the loading direction as shown in Figure 77. 81

In a similar analogy to demonstrate the non-directionality of the equi-biaxial fracture strain, the anisotropic ductile fracture criterion is proposed in the combination of the anisotropic stress triaxiality based on the Yld91 criterion as 

1 Y

C1 (θS ) *

ηA +

3 − LP + CO 6Y 1 + CO

+!C2 (θS )

!C4 (θS )

1 2 + ηA

5 9

ε¯pf = C3 (θS )

(A.26)

It can be verified that the proposed formula is able to consider the non-directionality of the equi-biaxial fracture strain from the fact that Y becomes unity particularly when the material is subjected to the equi-biaxial stress state under the plane stress condition. In addition, the model parameter C3 (θS ) in Equation A.26 still has the same meaning of ε¯bf as before. If the non-directionality of the equi-biaxial fracture strain is considered in the calibration of the model parameters by eliminating the directional dependence of the parameter C3 (θS ) through P3 = Q3 = R3 = ε¯bf . The final form of the uncoupled anisotropic ductile fracture criterion can be rewritten with the postulate of the plane stress condition for sheet metal forming application as follows: 

1 Y0

C1 (θS ) *

ηA +

3 − LP + CO 6Y 0 1 + CO

+!C2 (θS )

!C4 (θS )

1 2 + ηA

5 9

ε¯pf = C3 (θS )

(A.27)

It is worth to note that the proposed fracture criterion and its reduced one defined by the Hill’s 48 criterion exhibit the same performance in predicting the fracture strains used in the calibration of the model parameters due to its empirical characteristic resulting from the present form of each formula. It is, therefore, crucial to employ the measure of equivalent plastic strain which adequately demonstrates the deformation behavior of the tested material. It is also important to use its corresponding anisotropic stress triaxiality definition for the reliable fracture prediction. The generalized anisotropic stress triaxiality is capable of demonstrating the deformation behavior of various materials since the Yld91 criterion can reduce to various yield criteria including the Hill’s 48, Hosford, von Mises, and Tresca criteria. This implies that the proposed fracture criterion has a considerable potential to deal with fracture behavior of various materials in the same explicit form in addition to the characteristic of non-directionality of the equi-biaxial fracture strain. In general, there exists the constitutive law to describe the yielding of material, which has an influence on the shape of the fracture forming limit diagram. It is, thus, necessary to have an in-depth understanding about which yield criterion is suitable for adequately demonstrating the material behavior during the sheet metal forming process. After figuring out one of the most suitable constitutive laws for the targeted material, it should be also 82

applied to the fracture criterion by the anisotropic stress triaxiality defined by the constitutive law in consideration of the consistency in modeling of both fracture and yield criteria.

83

Appendix B. Construction of the Anisotropic Fracture Forming Limit Diagram Based on the Non-quadratic Yield Criterion A theoretical transformation process is necessary to construct the FFLD from the uncoupled ductile fracture criterion. In general, a transformation law is explicitly derived from the corresponding work-conjugate of the equivalent stress under the proportional loading. It is, however, not allowed to make an explicit form of the corresponding work-conjugate for nonquadratic yield criteria. As an alternative way to construct the FFLD, the explicit relation between the major and minor strains is rigorously derived in terms of the strain and stress paths based on the plastic work equivalence. It is worth to note that both methods to construct the FFLD are, in fact, equivalent from each other because those methods have the same base of the plastic work equivalence in their way of construction. The plastic work equivalence is defined as ˙ = σij ε˙p = σi ε˙p = σ1 ε˙p + σ2 ε˙p + σ3 ε˙p = σ W ¯A ε¯˙p 1 2 3 ij i

(B.0)

Assuming the plane stress condition (σ3 = 0) yields   p σ2 σ ¯A ε¯˙ ε˙p2 σ2 ε˙p2 p p ˙ = σ ¯ , β= ε ¯ → ε ˙ = where α = σ1 ε˙p1 1 + A 1 σ1 ε˙p1 σ1 1 + αβ σ1 ε˙p1

(B.1)

When the proportional loading condition is imposed to Equation B.1, the relation between the major strain and the minor strain is defined as Z 0

tf

ε˙p1 dt =

Z 0

tf

p p ˜ /∂σ2 p σ ¯A ε¯f ∂g (σ) σ ¯A ε¯˙ dt ⇒ εp1 = , εp2 = αεp1 = ε ˜ /∂σ1 1 σ1 1 + αβ σ1 1 + αβ ∂g (σ)

(B.2)

˜ represent time elapsed to fracture, the equivalent plastic strain at the onwhere tf , ε¯pf , and g(σ) set of fracture predicted from the fracture criterion and a plastic potential, respectively. An explicit expression of σ ¯A /σi (i = 1, 2, 3) can be obtained in terms of the anisotropic stress triaxiality and the Lode parameter by following the transformation methodology as discussed in Appendix A.1 for the non-quadratic yield criterion. If the plane stress is assumed (σi = 0, with i = 1, 2, 3) for sheet metal application, the relation between the anisotropic stress triaxiality and the Lode parameter is uniquely defined. This makes possible to explicitly express σ ¯A /σi (i = 1, 2, 3) and 84

the anisotropic ductile fracture criterion ε¯pf as a function of the Lode parameter only. Since the Lode parameter can be expressed by the ratio of minor to major stress under the plane stress condition, we can define the terms in Equation B.2 by the stress ratio except for α. Note that the strain ratio can be described by the stress ratio by associative flow rule. Following the transformation procedure introduced, an explicit relation between the major and minor strains is rigorously obtained in consideration of the general principal stress states as follows: Intermediate triaxiality (From the uniaxial tension to the equi-biaxial tension) (σ1 > 0, σ2 > 0, σ3 = 0) ε1 = ε¯pf

σ ¯A (1 + βα) = ε¯pf Y (1 + βα) , ε2 = αεp1 σ1

where β =

σ2 1 + LP = σ1 2

(B.3a)

Low and negative stress triaxiality (From the uniaxial compression to the uniaxial tension) (σ1 > 0, σ2 = 0, σ3 < 0) ε1 =

σ ¯A ε¯pf σ1

(1 + βα) =

ε¯pf



2Y 1 − LP



(1 + βα) , ε2 = αεp1

where β =

σ3 1 + LP =− σ1 1 − LP

(B.3b)

Negative stress triaxiality (From the equi-biaxial compression to the uniaxial compression) (σ1 = 0, σ2 < 0, σ3 < 0) ε1 =

σ ¯A ε¯pf σ2

(1 + βα) =

ε¯pf



2Y −1 + LP



(1 + βα) , ε2 = αεp1

where β =

σ3 2 = σ2 1 − LP

(B.3c)

where ε¯pf

 = C3 (θS )

1 Y

−C1 (θS ) *

ηA +

3 − LP + CO 6Y 1 + CO

+!−C2 (θS )

!−C4 (θS )

1 2 + ηA

5 9

(B.3d)

The explicit form of a function Y varies with the non-quadratic yield criterion employed in the construction of the anisotropic fracture criterion as follows: For the Yld2004-18p criterion 1 )   1 ( 1,3  0   00  0 − B 00 m m √ θ + 2π (i − 2) θ + 2π (j − 2) B 1 m P √ A 0 Y = − 2 −A00 cos A + 2 −A cos 4 3 3 3 i,j

(B.4a)

For the Yld91 criterion   1 −1      m     1 p θA θA + π m m 1 m π − θA m sin + sin + sin Y = −3Q 2 3 3 3

(B.4b)

In Equation B.3, the relation between the strain and stress paths is defined by the associative flow rule. Consequently, the anisotropic FFLD can be constructed according to the loading 85

direction when the non-quadratic yield criterion is employed. Note that the overall shapes of the anisotropic FFLD as well as the fracture locus are strongly dependent on the exponent m which determines the curvature of the Yld2004-18p and the Yld91 criteria. The change of yield surface curvature affects the value of anisotropic stress triaxiality, which eventually varies the shape of the fracture forming limit. For example, Figures 78 and 79 represent the influence of the exponent m on the yield surface and the fracture forming limits for the loading direction of 0°. Here, the parameters ci of the Yld91 criterion are set to the values in Table 29.

86

Bibliography [1] Report: Advanced high-strength steels application guidelines, Section 7 – Appendix case studies. World Auto Steel. [2] Standard test method for plastic strain ratio r for sheet metal. ASTM Standards E51700(2010), 2010. doi: 10.1520/E0517-00R10, www.astm.org. [3] Aramis v6.3.0 User manual software: Bulge test evalutation (biaxial yield stress calculation). GOM mbH, 2014. [4] H. Aretz. A non-quadratic plane stress yield function for orthotropic sheet metals. Journal of Materials Processing Technology, 168(1):1–9, 2005. [5] R. Arrioux, C. Bedrin, and M. Boivin. Determination of an intrinsic forming limit stress diagram for isotropic metal sheets. In International Deep Drawing Research Group, Proceedings of the 12th Biennial Congress IDDRG, pages 61–71, 1982. [6] M. Atkinson. Accurate determination of biaxial stress–strain relationships from hydraulic bulging tests of sheet metals. International Journal of Mechanical Sciences, 39(7):761– 769, 1997. [7] G. Bae. Strain rate-dependent flow stress curves in the large deformation range. In NUMISHEET 2014: The 9th International Conference and Workshop on Numerical Simulation of 3D Sheet Metal Forming Processes: Part A Benchmark Problems and Results and Part B General Papers, volume 1567, pages 504–507. AIP Publishing, 2013. [8] Y. Bai. Effect of loading history in necking and fracture. Ph.D. dissertation, Massachusetts Institute of Technology, 2008. [9] Y. Bai and T. Wierzbicki. Forming severity concept for predicting sheet necking under complex loading histories. International Journal of Mechanical Sciences, 50(6):1012–1022, 2008. [10] Y. Bai and T. Wierzbicki. Application of extended Mohr–Coulomb criterion to ductile fracture. International Journal of Fracture, 161(1):1–20, 2010.

87

[11] D. Banabic. Sheet metal forming processes: Constitutive modelling and numerical simulation. Springer, 2010. [12] D. Banabic. Advances in plastic anisotropy and forming limits in sheet metal forming. In ASME 2015 International Manufacturing Science and Engineering Conference, pages V001T02A090–V001T02A090. American Society of Mechanical Engineers, 2015. [13] D. Banabic, H. Aretz, D. S. Comsa, and L. Paraianu. An improved analytical description of orthotropy in metallic sheets. International Journal of Plasticity, 21(3):493–512, 2005. [14] D. Banabic, T. Balan, and D. S. Comsa. A new yield criterion for orthotropic sheet metals under plane-stress conditions. In Proceedings of the 7th Conference ‘TPR2000’, Cluj Napoca, Romania, pages 217–224, 2000. [15] D. Banabic, T. Kuwabara, T. Balan, D. S. Comsa, and D. Julean. Non-quadratic yield criterion for orthotropic sheet metals under plane-stress conditions. International Journal of Mechanical Sciences, 45(5):797–811, 2003. [16] D. Banabic, W. M¨ uller, and K. P¨ohlandt. Determination of yield loci from cross tensile tests assuming various kinds of yield criteria. Sheet metal forming beyond, pages 343–349, 2000. [17] Y. Bao and T. Wierzbicki. On fracture locus in the equivalent strain and stress triaxiality space. International Journal of Mechanical Sciences, 46(1):81–98, 2004. [18] F. Barlat, H. Aretz, J. Yoon, M. Karabin, J. Brem, and R. Dick. Linear transfomationbased anisotropic yield functions. International Journal of Plasticity, 21(5):1009–1039, 2005. [19] F. Barlat, H. Aretz, J. Yoon, M. Karabin, J. Brem, and R. Dick. Linear transfomationbased anisotropic yield functions. International Journal of Plasticity, 21(5):1009–1039, 2005. [20] F. Barlat, J. C. Brem, J. W. Yoon, K. Chung, R. E. Dick, D. J. Lege, F. Pourboghrat, S. H. Choi, and E. Chu. Plane stress yield function for aluminum alloy sheets – Part 1: Theory. International Journal of Plasticity, 19(9):1297–1319, 2003. [21] F. Barlat, D. J. Lege, and J. C. Brem. A six-component yield function for anisotropic materials. International Journal of Plasticity, 7(7):693–712, 1991.

88

[22] F. Barlat and K. Lian. Plastic behavior and stretchability of sheet metals. Part I: A yield function for orthotropic sheets under plane stress conditions. International Journal of Plasticity, 5(1):51–66, 1989. [23] F. Barlat, Y. Maeda, K. Chung, M. Yanagawa, J. C. Brem, Y. Hayashida, D. J. Lege, K. Matsui, S. J. Murtha, S. Hattori, R. C. Becker, and S. Makosey. Yield function development for aluminum alloy sheets. Journal of the Mechanics and Physics of Solids, 45(11-12):1727–1763, 1997. [24] I. Barsoum and J. Faleskog.

Rupture mechanisms in combined tension and shear–

Experiments. International Journal of Solids and Structures, 44(6):1768–1786, 2007. [25] I. Barsoum and J. Faleskog.

Rupture mechanisms in combined tension and shear–

Micromechanics. International Journal of Solids and Structures, 44(17):5481–5498, 2007. [26] K.-J. Bathe. Finite element procedures. Klaus-Jurgen Bathe, 2006. [27] A. M. Beese, M. Luo, Y. Li, Y. Bai, and T. Wierzbicki. Partially coupled anisotropic fracture model for aluminum sheets. Engineering Fracture Mechanics, 77(7):1128–1152, 2010. [28] J. Bishop and R. Hill. XLVI. A theory of the plastic distortion of a polycrystalline aggregate under combined stresses. The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science, 42(327):414–427, 1951. [29] T. Børvik, O. S. Hopperstad, and K. O. Pedersen. Quasi-brittle fracture during structural impact of AA7075-T651 aluminium plates. International Journal of Impact Engineering, 37(5):537–551, 2010. [30] P. Bridgman. The effect of hydrostatic pressure on the fracture of brittle substances. Journal of Applied Physics, 18(2):246–258, 1947. [31] P. Bridgman. Studies in large plastic flow and fracture with special emphasis on the effects of hydrostatic pressure. 1st ed. McGraw-Hill, New York, 1952. [32] F. Bron and J. Besson. A yield function for anisotropic materials application to aluminum alloys. International Journal of Plasticity, 20(4):937–963, 2004. [33] P. Brozzo, B. Deluca, and R. Rendina. A new method for the prediction of formability limits in metal sheets. In Proc. 7th biennal Conf. IDDR, 1972. 89

[34] M. Br¨ unig. An anisotropic ductile damage model based on irreversible thermodynamics. International Journal of Plasticity, 19(10):1679–1713, 2003. [35] M. Br¨ unig. Numerical analysis of anisotropic ductile continuum damage. Computer Methods in Applied Mechanics and Engineering, 192(26):2949–2976, 2003. [36] M. Br¨ unig. Continuum framework for the rate-dependent behavior of anisotropically damaged ductile metals. Acta mechanica, 186(1):37–53, 2006. [37] M. Br¨ unig, D. Albrecht, and S. Gerke. Modeling of ductile damage and fracture behavior based on different micromechanisms. International Journal of Damage Mechanics, 20(4):558–577, 2011. [38] B. Budiansky, J. Hutchinson, and S. Slutsky. Void growth and collapse in viscous solids. Pergamon Press, London, 1982. [39] R. P. Cardoso and J. W. Yoon.

Stress integration method for a nonlinear kine-

matic/isotropic hardening model and its characterization based on polycrystal plasticity. International Journal of Plasticity, 25(9):1684–1710, 2009. [40] O. Cazacu and F. Barlat. Generalization of Drucker’s yield criterion to orthotropy. Mathematics and Mechanics of Solids, 6(6):613–630, 2001. [41] O. Cazacu and F. Barlat. Application of the theory of representation to describe yielding of anisotropic aluminum alloys. International Journal of Engineering Science, 41(12):1367– 1385, 2003. [42] O. Cazacu, B. Plunkett, and F. Barlat. Orthotropic yield criterion for hexagonal closed packed metals. International Journal of Plasticity, 22(7):1171–1194, 2006. [43] J. L. Chaboche. Continuous damage mechanics: A tool to describe phenomena before crack initiation. Nuclear Engineering and Design, 64(2):233–247, 1981. [44] J. L. Chaboche. Continuum damage mechanics: Part I – general concepts. Journal of Applied Mechanics, 55(1):59–64, 1988. [45] J. L. Chaboche. Continuum damage mechanics: Part II – damage growth, crack initiation, and crack growth. Journal of Applied Mechanics, 55(1):65–72, 1988. [46] W. F. Chen and D. J. Han. Plasticity for structural engineers. J. Ross Publishing, 2007.

90

[47] F. Chinesta, E. Cueto, D. Banabic, F. Barlat, O. Cazacu, and T. Kuwabara. Anisotropy and formability. In Advances in Material Forming, pages 143–173. Springer Paris, 2007. [48] T. C. Chu, W. F. Ranson, and M. A. Sutton. Applications of digital image correlation techniques to experimental mechanics. Experimental Mechanics, 25(3):232–244, 1985. [49] S. E. Clift, P. Hartley, C. Sturgess, and G. Rowe. Fracture prediction in plastic deformation processes. International Journal of Mechanical Sciences, 32(1):1–17, 1990. [50] M. Cockcroft and D. Latham. Ductility and the workability of metals. Journal of the Institute of Metals, 96(1):33–39, 1968. [51] D. S. Comsa and D. Banabic. Numerical simulation of sheet metal forming processes using a new yield criterion. Key Engineering Materials, 344:833–840, 2007. [52] D. S. Comsa and D. Banabic. Plane stress yield criterion for highly anisotropic sheet metals. In Numisheet 2008, Interlaken, Switzerland, pages 43–48, 2008. [53] L. Cortese, F. Nalli, and M. Rossi. A nonlinear model for ductile damage accumulation under multiaxial non-proportional loading conditions. International Journal of Plasticity, 85:77–92, 2016. [54] M. Crisfield and V. Ciampi. Non-linear finite element analysis of solids and structures vol. 1. Meccanica, 32(6):586–587, 1997. ˇ Lozina. A finite element formulation based on non-associated [55] V. Cvitani´c, F. Vlak, and Z. plasticity for sheet metal forming. International Journal of Plasticity, 24(4):646–687, 2008. [56] K. Danas and P. P. Casta˜ neda. Influence of the lode parameter and the stress triaxiality on the failure of elasto-plastic porous materials. International Journal of Solids and Structures, 49(11):1325–1342, 2012. [57] R. O. Davis and A. P. Selvadurai. Plasticity and geomechanics. Cambridge university press, 2005. [58] D. C. Drucker. Some implications of work hardening and ideal plasticity. Quarterly of Applied Mathematics, 7(4):411–418, 1950. [59] D. C. Drucker. A more fundamental approach to plastic stress-strain relations. In Journal of Applied Mechanics-Transactions of the ASME, volume 18, pages 323–323. ASMEAMER Soc Mechanical Eng 345 E 47th St, NEW YORK, NY 10017, 1951. 91

[60] D. C. Drucker. A definition of stable inelastic material. Brown Univ Providence RI, 1957. [61] D. C. Drucker and W. Prager. Soil mechanics and plastic analysis or limit design. Quarterly of applied mathematics, 10(2):157–165, 1952. [62] M. Dunand, A. P. Maertens, M. Luo, and D. Mohr. Experiments and modeling of anisotropic aluminum extrusions under multi-axial loading – Part I: Plasticity. International Journal of Plasticity, 36(0):34–49, 2012. [63] M. Dunand and D. Mohr. Hybrid experimental–numerical analysis of basic ductile fracture experiments for sheet metals. International Journal of Solids and Structures, 47(9):1130– 1143, 2010. [64] C. Efstathiou, H. Sehitoglu, and J. Lambros. Multiscale strain measurements of plastically deforming polycrystalline titanium: Role of deformation heterogeneities. International Journal of Plasticity, 26(1):93–106, 2010. [65] A. Eringen. Tensor analysis. Continuum physics, 1:1–155, 1971. [66] A. C. Eringen. Nonlinear theory of continuous media. McGraw-Hill, 1962. [67] G. Ferron, R. Makkouk, and J. Morreale. A parametric description of orthotropic plasticity in metal sheets. International Journal of Plasticity, 10(5):431–449, 1994. [68] A. Freudenthal. The inelastic behavior of solids. Wiley, New York, 1950. [69] X. Gao, T. Wang, and J. Kim. On ductile fracture initiation toughness: Effects of void volume fraction, void shape and void distribution. International Journal of Solids and Structures, 42(18):5097–5117, 2005. [70] X. Gao, G. Zhang, and C. Roe. A study on the effect of the stress state on ductile fracture. International Journal of Damage Mechanics, 2009. [71] H. Ghadbeigi, C. Pinna, S. Celotto, and J. Yates. Local plastic strain evolution in a high strength dual-phase steel. Materials Science and Engineering: A, 527(18–19):5026–5032, 2010. [72] A. Godara and D. Raabe. Influence of fiber orientation on global mechanical behavior and mesoscale strain localization in a short glass-fiber-reinforced epoxy polymer composite during tensile deformation investigated using digital image correlation. Composites science and technology, 67(11):2417–2427, 2007. 92

[73] M. Gologanu, J. B. Leblond, and J. Devaux. Approximate models for ductile metals containing non-spherical voids–case of axisymmetric prolate ellipsoidal cavities. Journal of the Mechanics and Physics of Solids, 41(11):1723–1754, 1993. [74] G. M. Goodwin. Application of strain analysis to sheet metal forming problems in the press shop. SAE technical paper, 1968. [75] M. Gotoh. A theory of plastic anisotropy based on a yield function of fourth order (plane stress state)–I. International Journal of Mechanical Sciences, 19(9):505–512, 1977. [76] K.-H. Grote and E. K. Antonsson. Springer handbook of mechanical engineering. Springer Science & Business Media, 10, 2009. [77] G. Gu and D. Mohr. Anisotropic Hosford–Coulomb fracture initiation model: Theory and application. Engineering Fracture Mechanics, 147:480–497, 2015. [78] A. L. Gurson et al. Continuum theory of ductile rupture by void nucleation and growth: Part I – yield criteria and flow rules for porous ductile media. Journal of engineering materials and technology, 99(1):2–15, 1977. ˇ [79] M. Haliloviˇc, M. Vrh, and B. Stok. Nice – An explicit numerical scheme for efficient integration of nonlinear constitutive equations. Mathematics and Computers in Simulation, 80(2):294–313, 2009. [80] H. Hayashi and T. Nakagawa. Recent trends in sheet metals and their formability in manufacturing automotive panels. Journal of Materials Processing Technology, 46(3):455– 487, 1994. [81] S. Hecker. Experimental studies of yield phenomena in biaxially loaded metals. Constitutive equations in viscoplasticity: computational and engineering aspects, pages 1–33, 1976. [82] A. Hershey. The plasticity of an isotropic aggregate of anisotropic face centered cubic crystals. Journal of Applied Mechanics, 21:241–249, 1954. [83] R. Hill. A theory of the yielding and plastic flow of anisotropic metals. Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 193(1033):281– 297, 1948. [84] R. Hill. Theoretical plasticity of textured aggregates. Mathematical Proceedings of the Cambridge Philosophical Society, 85(01):179–191, 1979. 93

[85] R. Hill. Constitutive modeling of orthotropic plasticity in sheet metals. Journal of the Mechanics and Physics of Solids, 38(3):405–417, 1990. [86] R. Hill. A user-friendly theory of orthotropic plasticity in sheet metals. International Journal of Mechanical Sciences, 35(1):19–25, 1993. [87] R. t. Hill. On discontinuous plastic states, with special reference to localized necking in thin sheets. Journal of the Mechanics and Physics of Solids, 1(1):19–30, 1952. [88] H. Hippke, N. Manopulo, J. W. Yoon, and P. Hora. On the efficiency and accuracy of stress integration algorithms for constitutive models based on non-associated flow rule. International Journal of Material Forming, pages 1–8, 2017. [89] P. Hora, L. Tong, and J. Reissner. A prediction method for ductile sheet metal failure in fe-simulation. In Proceedings of NUMISHEET, volume 96, pages 252–256, 1996. [90] W. Hosford. On yield loci of anisotropic cubic metals. In Proc. 7th North American Metalworking Conf., SME, Dearborn, MI, pages 191–197, 1979. [91] W. F. Hosford. A generalized isotropic yield criterion. Journal of Applied Mechanics, 39(2):607–609, 1972. [92] W. Hu. An orthotropic yield criterion in a 3D general stress state. International Journal of Plasticity, 21(9):1771–1796, 2005. [93] J. Huh, H. Huh, and C. S. Lee. Effect of strain rate on plastic anisotropy of advanced high strength steel sheets. International Journal of Plasticity, 44:23–46, 2013. [94] L. Kachanov. Time of the rupture process under creep conditions, izu. Akad. Nauk SSR Otd. Tech, pages 26–31, 1958. [95] A. P. Karafillis and M. C. Boyce. A general anisotropic yield criterion using bounds and a transformation weighting tensor. Journal of Materials Processing Technology, 41:1859– 1886, 1993. [96] S. P. Keeler and W. A. Backofen. Plastic instability and fracture in sheets stretched over rigid punches. ASM Transactions Quarterly, 56(1):25–48, 1963. [97] A. S. Khan and S. Huang. Continuum theory of plasticity. John Wiley & Sons, 1995. [98] A. S. Khan and H. Liu. A new approach for ductile fracture prediction on Al2024–T351 alloy. International Journal of Plasticity, 35:1–12, 2012. 94

[99] A. S. Khan and H. Liu. Strain rate and temperature dependent fracture criteria for isotropic and anisotropic metals. International Journal of Plasticity, 37:1–15, 2012. [100] J. Kim, X. Gao, and T. Srivatsan. Modeling of crack growth in ductile solids: A threedimensional analysis. International journal of solids and structures, 40(26):7357–7374, 2003. [101] J. Kim, X. Gao, and T. S. Srivatsan. Modeling of void growth in ductile solids: Effects of stress triaxiality and initial porosity. Engineering Fracture Mechanics, 71(3):379–400, 2004. [102] M. S. Kirugulige, H. V. Tippur, and T. S. Denney. Measurement of transient deformations using digital image correlation method and high-speed photography: Application to dynamic fracture. Applied Optics, 46(22):5083–5096, 2007. [103] Y. K. Ko, J. S. Lee, H. Huh, H. Kim, and S. Park. Prediction of fracture in hub-hole expanding process using a new ductile fracture criterion. Journal of materials processing technology, 187:358–362, 2007. ¨ N. Cora. An experimental study on the comparative assessment [104] M. Ko¸c, E. Billur, and O. of hydraulic bulge test analysis methods. Materials & Design, 32(1):272–281, 2011. [105] Y. P. Korkolis and S. Kyriakides. Inflation and burst of aluminum tubes. Part II: An advanced yield function including deformation-induced anisotropy. International Journal of Plasticity, 24(9):1625–1637, 2008. [106] P. V. Lade, R. B. Nelson, and Y. M. Ito. Non-associated flow and stability of granular materials. Journal of Engineering Mechanics, 113(9):1302–1318, 1987. [107] W. Lankford, S. Snyder, and J. Bauscher. New criteria for predicting the press performance of deep drawing sheets. Transactions of the American Society of Metals, 42:1197– 1232, 1950. [108] A. G. Leacock. A mathematical description of orthotropy in sheet metals. Journal of the Mechanics and Physics of Solids, 54(2):425–444, 2006. [109] J. Lemaitre. A continuous damage mechanics model for ductile fracture. Transactions of the ASME. Journal of Engineering Materials and Technology, 107(1):83–89, 1985. [110] J. Lemaitre. A course on damage mechanics. Springer Science & Business Media, 2012. 95

[111] J. Lemaitre, R. Desmorat, and M. Sauzay. Anisotropic damage law of evolution. European Journal of Mechanics-A/Solids, 19(2):187–208, 2000. [112] D. Li, Z. Chen, L. Sun, J. Lee, and R. Wagoner. An improved test for shear fracture. International Journal of Solids and Structures, 97:29–42, 2016. [113] Y. Li, M. Luo, J. Gerlach, and T. Wierzbicki. Prediction of shear-induced fracture in sheet metal forming. Journal of Materials Processing Technology, 210(14):1858–1869, 2010. [114] S. Lin and J. Ding. A modified form of Hill’s orientation-dependent yield criterion for orthotropic sheet metals. Journal of the Mechanics and Physics of Solids, 44(11):1739– 1764, 1996. [115] A. Lipp, I. Heinle, P. Craighero, H. Friebe, M. Klein, and T. M¨oller. Improved flow curve determination using the bulge test combined with optical measurement systems and compensation strategies. In IDDRG, volume 2012, pages 224–231, 2012. [116] W. Lode. Versuche u ¨ber den einfluß der mittleren hauptspannung auf das fließen der metalle eisen, kupfer und nickel. Zeitschrift f¨ ur Physik, 36(11-12):913–939, 1926. [117] R. W. Logan and W. F. Hosford. Upper-bound anisotropic yield locus calculations assuming -pencil glide. International Journal of Mechanical Sciences, 22(7):419–430, 1980. [118] Y. Lou and H. Huh. Extension of a shear-controlled ductile fracture model considering the stress triaxiality and the lode parameter. International Journal of Solids and Structures, 50(2):447–455, 2013. [119] Y. Lou, H. Huh, S. Lim, and K. Pack. New ductile fracture criterion for prediction of fracture forming limit diagrams of sheet metals. International Journal of Solids and Structures, 49(25):3605–3615, 2012. [120] Y. Lou, J. W. Yoon, and H. Huh. Modeling of shear ductile fracture considering a changeable cut-off value for stress triaxiality. International Journal of Plasticity, 54:56– 80, 2014. [121] J. Lubliner. Plasticity theory. Courier Corporation, 2008. [122] A. C. Luo. Constitutive laws and damage theory. In Nonlinear Deformable-body Dynamics, pages 135–160. Springer, 2010. 96

[123] M. Luo, M. Dunand, and D. Mohr. Experiments and modeling of anisotropic aluminum extrusions under multi-axial loading – Part II: Ductile fracture. International Journal of Plasticity, 32–33(0):36–58, 2012. [124] M. Luo and T. Wierzbicki. Numerical failure analysis of a stretch-bending test on dualphase steel sheets using a phenomenological fracture model. International Journal of Solids and Structures, 47(22):3084–3102, 2010. [125] Z. Marciniak and K. Kuczy´ nski. Limit strains in the processes of stretch-forming sheet metal. International journal of mechanical sciences, 9(9):609IN1613–612IN2620, 1967. [126] J. E. Marsden and T. J. Hughes. Mathematical foundations of elasticity. Courier Corporation, 1994. [127] P. Martins, N. Bay, A. Tekkaya, and A. Atkins. Characterization of fracture loci in metal forming. International Journal of Mechanical Sciences, 83:112–123, 2014. [128] F. A. McClintock. A criterion for ductile fracture by the growth of holes. Journal of applied mechanics, 35(2):363–371, 1968. [129] J. Min, J. E. Carsley, J. Lin, Y. Wen, and B. Kuhlenk¨otter. A non-quadratic constitutive model under non-associated flow rule of sheet metals with anisotropic hardening: Modeling and experimental validation. International Journal of Mechanical Sciences, 119:343–359, 2016. [130] J. Min, L. G. Hector Jr, J. E. Carsley, T. B. Stoughton, B. E. Carlson, and J. Lin. Spatiotemporal characteristics of plastic instability in AA5182-O during biaxial deformation. Materials & Design, 83:786–794, 2015. [131] J. Min, T. B. Stoughton, J. E. Carsley, B. E. Carlson, J. Lin, and X. Gao. Accurate characterization of biaxial stress–strain response of sheet metal from bulge testing. International Journal of Plasticity, 2016. [132] D. Mohr and S. J. Marcadet. Micromechanically-motivated phenomenological Hosford– Coulomb model for predicting ductile fracture initiation at low stress triaxialities. International Journal of Solids and Structures, 67:40–55, 2015. [133] K. Nahshon and J. Hutchinson. Modification of the Gurson model for shear failure. European Journal of Mechanics-A/Solids, 27(1):1–17, 2008.

97

[134] A. Nasser, A. Yadav, P. Pathak, and T. Altan. Determination of the flow stress of five AHSS sheet materials (DP 600, DP 780, DP 780-CR, DP 780-HY and TRIP 780) using the uniaxial tensile and the biaxial viscous pressure bulge (VPB) tests. Journal of Materials Processing Technology, 210(3):429–436, 2010. [135] K. L. Nielsen and V. Tvergaard. Ductile shear failure or plug failure of spot welds modelled by modified Gurson model. Engineering Fracture Mechanics, 77(7):1031–1047, 2010. [136] S. Oh, C. Chen, and S. Kobayashi. Ductile fracture in axisymmetric extrusion and drawing – Part 2: Workability in extrusion and drawing. Journal of Engineering for Industry, 101(1):36–44, 1979. [137] N. S. Ottosen and M. Ristinmaa. The mechanics of constitutive modeling. Elsevier, 2005. [138] M. Oyane, T. Sato, K. Okimoto, and S. Shima. Criteria for ductile fracture and their applications. Journal of Mechanical Working Technology, 4(1):65–81, 1980. [139] K. Pack and C. C. Roth. The second Sandia Fracture Challenge: Blind prediction of dynamic shear localization and full fracture characterization. International Journal of Fracture, 198(1-2):197–220, 2016. [140] J. Papasidero, V. Doquet, and D. Mohr. Ductile fracture of aluminum 2024-T351 under proportional and non-proportional multi-axial loading: Bao–Wierzbicki results revisited. International Journal of Solids and Structures, 69:459–474, 2015. [141] L. Paraianu, D. S. Comsa, G. Cosovici, P. Jurco, and D. Banabic. An improvement of the BBC2000 yield criterion. In Proceedings of the ESAFORM 2003 Conference, Salerno, 2002. [142] T. Park and K. Chung. Non-associated flow rule with symmetric stiffness modulus for isotropic-kinematic hardening and its application for earing in circular cup drawing. International Journal of Solids and Structures, 49(25):3582–3593, 2012. [143] E. Parsons, M. C. Boyce, and D. M. Parks. An experimental investigation of the largestrain tensile behavior of neat and rubber-toughened polycarbonate. Polymer, 45(8):2665– 2684, 2004. [144] R. Pearce. Some aspects of anisotropic plasticity in sheet metals. International Journal of Mechanical Sciences, 10(12):995–1004, 1968.

98

[145] W. H. Peters, W. F. Ranson, M. A. Sutton, T. C. Chu, and J. Anderson. Applications of digital image correlation methods to rigid body mechanics. Optical Mechanics, 22:738– 742, 1983. [146] P. Prates, M. Oliveira, and J. Fernandes. Identification of material parameters for thin sheets from single biaxial tensile test using a sequential inverse identification strategy. International Journal of Material Forming, 9(4):547–571, 2016. [147] L.-Y. Qian, G. Fang, and P. Zeng. Modeling of the ductile fracture during the sheet forming of aluminum alloy considering non-associated constitutive characteristic. International Journal of Mechanical Sciences, 126:55–66, 2017. [148] W. F. Ranson and W. H. Peters. Digital imaging techniques in experimental stress analysis. Optical Engineering, 21(3):427–431, 1982. [149] A. Ranta-Eskola. Use of the hydraulic bulge test in biaxial tensile testing. International Journal of Mechanical Sciences, 21(8):457–465, 1979. [150] J. R. Rice and D. M. Tracey. On the ductile enlargement of voids in triaxial stress fields. Journal of the Mechanics and Physics of Solids, 17(3):201–217, 1969. [151] O. Rodriguez. Des lois geometriques qui regissent les desplacements d’un systeme solide dans l’espace et de la variation des coordonnees provenant de deplacements consideres independamment des causes qui peuvent les produire. Journal de Math´ematiques Pures et Appliqu´ees, 5:380–440, 1840. [152] G. Rousselier, M. Luo, and D. Mohr. Macroscopic plasticity modeling of anisotropic aluminum extrusions using a reduced texture methodology. International Journal of Plasticity, 30–31:144–165, 2012. [153] K. Saanouni and J. Chaboche. Computational damage mechanics. application to metal forming. Numerical and Computational Methods in Comprehensive Structural Integrity, 3:321–376, 2003. ˇ [154] F. Sebek. Ductile fracture criteria in multiaxial loading – Theory, experiments and application. PhD dissertation, Brno University of Technology, 2016. [155] H.-C. Shih. An experimental study on shear fracture of advanced high strength steels. United States Steel Corporation.

99

[156] J. C. Simo and T. J. Hughes. Computational inelasticity. Springer Science & Business Media, 7, 2006. [157] S. Soare, J. W. Yoon, and O. Cazacu. On the use of homogeneous polynomials to develop anisotropic yield functions with applications to sheet forming. International Journal of Plasticity, 24(6):915–944, 2008. [158] W. Spitzig and O. Richmond. The effect of pressure on the flow stress of metals. Acta metallurgica, 32(3):457–463, 1984. [159] R. V. Stone, T. Cox, J. Low, and J. Psioda. Microstructural aspects of fracture by dimpled rupture. International Metals Reviews, 30(1):157–180, 1985. [160] S. St¨oren and J. Rice. Localized necking in thin sheets. Journal of the Mechanics and Physics of Solids, 23(6):421–441, 1975. [161] T. B. Stoughton. A general forming limit criterion for sheet metal forming. International Journal of Mechanical Sciences, 42(1):1–27, 2000. [162] T. B. Stoughton. A non-associated flow rule for sheet metal forming. International Journal of Plasticity, 18(5):687–714, 2002. [163] T. B. Stoughton and J. W. Yoon. Sheet metal formability analysis for anisotropic materials under non-proportional loading. International journal of mechanical sciences, 47(12):1972–2002, 2005. [164] T. B. Stoughton and J. W. Yoon. Review of Drucker’s postulate and the issue of plastic stability in metal forming. International journal of plasticity, 22(3):391–433, 2006. [165] T. B. Stoughton and J. W. Yoon. Anisotropic hardening and non-associated flow in proportional loading of sheet metals. International Journal of Plasticity, 25(9):1777–1817, 2009. [166] T. B. Stoughton and J. W. Yoon. A new approach for failure criterion for sheet metals. International Journal of Plasticity, 27(3):440–459, 2011. [167] T. B. Stoughton and J. W. Yoon. Path-independent forming limits in strain and stress spaces. International Journal of Solids and Structures, 49(25):3616–3625, 2012. [168] T. B. Stoughton and X. Zhu. Review of theoretical models of the strain-based FLD and their relevance to the stress-based FLD. International Journal of Plasticity, 20(8):1463– 1486, 2004. 100

[169] H. Swift. Plastic instability under plane stress. Journal of the Mechanics and Physics of Solids, 1(1):1–18, 1952. [170] A. Taherizadeh, D. E. Green, A. Ghaei, and J.-W. Yoon. A non-associated constitutive model with mixed iso-kinematic hardening for finite element simulation of sheet metal forming. International journal of plasticity, 26(2):288–309, 2010. [171] V. Tarigopula, O. S. Hopperstad, M. Langseth, A. H. Clausen, F. Hild, O. G. Lademo, and M. Eriksson. A study of large plastic deformations in dual phase steel using digital image correlation and FE analysis. Experimental Mechanics, 48(2):197–197, 2008. [172] P. Thomason. Ductile fracture of metals. Pergamon Press plc, Ductile Fracture of Metals(UK), 1990,, page 219, 1990. [173] H. E. Tresca. M´emoire sur l’´ecoulement des corps solides. In M´emoires Present´es par Divers Savants, volume 18, pages 733–799. Acad´emie des Sciences, Paris, 1868. [174] V. Tvergaard and A. Needleman. Analysis of the cup-cone fracture in a round tensile bar. Acta metallurgica, 32(1):157–169, 1984. [175] H. Vegter, P. Drent, and J. Hu´etink. A planar isotropic yield criterion based on mechanical testing at multi-axial stress states. Simulation of Materials Processing: Theory, Methods and Applications, pages 345–350, 1995. [176] E. Voce. The relationship between stress and strain for homogeneous deformation. J Inst Met, 74:537–562, 1948. [177] R. von Mises. Mechanik der plastischen formanderung von kristallen. Zeitschrift fur Angewandte Mathematik und Mechanik, 8(3):161–185, 1928. [178] M. Vucetic, A. Bouguecha, I. Peshekhodov, T. G¨otze, T. Huinink, H. Friebe, T. M¨oller, and B.-A. Behrens. Numerical validation of analytical biaxial true stress. In The 8th International Conference and Workshop on numerical simulation of 3D seet metal forming processes (NUMISHEET 2011), pages 107–114, 2011. [179] R. H. Wagoner and G. R. Smith. Advanced high strength steel workshop. Arlington, Virginia, USA, pages 1–122, 2006. [180] T. Wierzbicki and L. Xue. On the effect of the third invariant of the stress deviator on ductile fracture. Impact and Crashworthiness Laboratory, Technical Report, (136), 2005. 101

[181] J. Woodthorpe and R. Pearce. The anomalous behaviour of aluminium sheet under balanced biaxial tension. International Journal of Mechanical Sciences, 12(4):341–347, 1970. [182] L. Xue. Damage accumulation and fracture initiation in uncracked ductile solids subject to triaxial loading. International Journal of Solids and Structures, 44(16):5163–5181, 2007. [183] L. Xue. Ductile fracture modeling – Theory, experimental investigation and numerical verification. PhD dissertation, Massachusetts Institute of Technology, 2007. [184] L. Xue. Constitutive modeling of void shearing effect in ductile fracture of porous materials. Engineering Fracture Mechanics, 75(11):3343–3366, 2008. [185] L. Xue and T. Wierzbicki. Ductile fracture initiation and propagation modeling using damage plasticity theory. Engineering Fracture Mechanics, 75(11):3276–3293, 2008. [186] J. Yoon, T. Stoughton, R. Dick, J. M. Cesar de Sa, and A. D. Santos. Earing prediction in cup drawing based on non-associated flow rule. In AIP Conference Proceedings, volume 908, pages 685–690. AIP, 2007. [187] J. W. Yoon, F. Barlat, R. E. Dick, K. Chung, and T. J. Kang. Plane stress yield function for aluminum alloy sheets – Part II: FE formulation and its implementation. International Journal of Plasticity, 20(3):495–522, 2004. [188] M. Yu. New system of strength theory. Xi’an Jiaotong University Press, Xi’an (in Chinese), 1992. [189] M. Yu and F. Liu. Smooth ridge model of generalized twin shear stress criterion [J]. Acta Mechanica Sinica, 2:012, 1990. [190] Z. Yue. Ductile damage prediction in sheet metal forming processes. PhD dissertation, Troyes, 2014. [191] D. Zeng, L. Chappuis, Z. C. Xia, and X. Zhu. A path independent forming limit criterion for sheet metal forming simulations. SAE International Journal of Materials & Manufacturing, 1(1):809–817, 2009. [192] G. Zhong-heng. Unified theory of variation principles in non-linear theory of elasticity. Applied Mathematics and Mechanics, 1(1):1–22, 1980. [193] X. Zhu, K. Weinmann, and A. Chandra. A unified bifurcation analysis of sheet metal forming limits. Transactions American Society of Mechanical Engineers Journal of Engineering Materials and Technology, 123(3):329–333, 2001. 102

[194] Y. Zhu, B. Dodd, R. Caddell, and W. Hosford. Convexity restrictions on non-quadratic anisotropic yield criteria. International Journal of Mechanical Sciences, 29(10-11):733– 741, 1987.

103

ε¯fp



2τmax σ ¯

d¯ εp  ε¯fpr η, θ¯

0

Z

ε¯fp

1 Cl3

Z 0

Cl1 

h1 + 3ηi 2

Cl2

d¯ εp



0 when x < 0

x when x ≥ 0

i 1h a a a (f1 − f2 ) + (f2 − f3 ) + (f1 − f3 ) 2

where hxi =

 1/n where ε¯fpr η, θ¯ = b(1 + c)

!−1/n

+ c (2η + f1 + f3 )

hπ hπ hπ i i i 2 2 2 cos 1 − θ¯ , f2 = cos 3 + θ¯ , f3 = − cos 1 + θ¯ 3 6 3 6 3 6

1/a

Z p ε ¯ f 1 σmax d¯ εp Ccl 0   Z p ε ¯ f 1 3σ max d¯ εp exp Crt 0 2¯ σ  Z ε¯fp  1 2 σmax d¯ εp Cb 0 3 σmax − σm Z ε¯fp 1 σ max d¯ εp Co 0 σ ¯ Z ε¯f   1 σm + Cos2 d¯ εp Cos1 0 σ ¯ Z ε¯fp 1 σ ¯ d¯ εp Cc 0 ( Z ε¯fp D x when x ≥ 0 1 σ σ E max m 1+3 d¯ εp where hxi = Ck 0 σ ¯ σ ¯ 0 when x < 0 Z ε¯fp d¯ εp  ε¯pr η, θ¯ 0 f ( " √  ¯  # "r ¯    ¯ #)−1/n  3 1 + c2 A θπ θπ 1 θπ 1 √ (cθax − cθs ) sec where ε¯fpr η, θ¯ = cθs + −1 cos + c1 η + sin c 6 3 6 3 6 2 − 3 2 (

Formulation of damage model

Table 1 Uncoupled isotropic ductile fracture criteria Remark

Cockcroft and Latham [50] Rice and Tracey [150] Brozzo et al. [33] Oh et al. [136] Oyane and Sato [138] Clift [49]

Ko et al. [103]

Bai and Wierzbicki [10]

Lou et al. [119]

Marcadet and Mohr [132]

f1 =

104

Formulation of damage model

Table 2 Uncoupled anisotropic ductile fracture criteria Remark

Luo et al. [123]

Gu and Mohr [77]

0

0

0

0

0

M22

0

0

0

1

0

0

0

M44

0

0

0

1

0

0

0 0 0

0

1

hπ hπ i i 2 2 cos 3 + θ¯ , f3 = − cos 1 + θ¯ 3 6 3 6  0  0    0    0   0 

r Z ε¯fp Z ε¯fp 2 β (dεp ) · β (dεp ) w (σ) d¯ εp = w (σ) 3 0 0 ( " √  ¯  # "r ¯    ¯ #)−1/n 2 θπ θπ 1 θπ 3 1 + c 1 √ (1 − c3 ) sec where w (σ) = c2 c3 + −1 cos + η + sin 6 3 6 3 6 2− 3   1 0 0 0 0 0    0 β22 0 0 0 0      0 0 0   0 0 β33 β=    0 β12 0 0   0 0   0 0 1 0   0 0 0 0 0 0 0 1 Z ε¯a,f  1/n d¯ ε 1 + c a pr ˆ σ] = b = 1 with ε¯a,f [σ/¯ pr ˆ σ] ε¯a,f [σ/¯ ˆ σ] gˆHC [Mσ/¯ 0  i1/a 1h a a a ˆ σ] = where gˆHC [σ/¯ + c (2η + f1 + f3 ) (f1 − f2 ) + (f2 − f3 ) + (f1 − f3 ) 2 hπ i 2 f1 = cos 1 − θ¯ , f2 = 3 6  1 M12 0 0 0      M=    

0

105

Table 3 Material quantities needed for the identification of several yield criteria [12] Author, year

σ0

σ30 σ45 σ75 σ90 σb

r0

r30 r45 r75 r90

rb

3D A1 A2

Hill’s family Hill 1948

x

Hill 1979

x

Hill 1990

x

Hill 1993

x x x

x

x

x

x

x

x

x

x

x

x

Lin, Ding 1996

x

x

x

x

Hu 2005

x

x

x

Leacock 2006

x

x

x

x x

x x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

Hershey’s family Hosford 1979

x

x

x

Barlat 1989

x

x

x

Barlat 1991

x

x

x

Karafillis, Boyce 1993

x

x

x

Barlat 1997

x

x

x

Banabic 2000

x

x

Barlat 2000

x

Bron, Besson 2003

x

Balat 2004

x

Banabic 2005

x

Banabic 2008

x

x

x

Cazacu–Barlat 2001

x

x

Cazacu–Barlat 2003

x

Cazacu–Pluncket 2006

x

x x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

x

Gotoh 1997

x

x

x

x

x

x

x

x

x

x

x

x

x

x

Comsa, Banabic 2007

x

x

x

x

x

x

x

x

x

x

Soare 2007 (Poly 4,6,8)

x

x

x

x

x

x

x

x

x

x

x

x

x

x

Drucker’s family

Polynomial criteria x x

x

x

x

x

x

Other criteria Ferron 1994

x

Vegter 1995

x

x x

x

x x

x

x x

x

x x

x

x

x

x

x

x

3D: the model is extendable to spatial stress states A1: the yield criterion is able to describe ‘the first order anomalous behavior’ A2: the yield criterion is able to describe ‘the second order anomalous behavior’

106

x

Table 4 Chemical compositiona of the DP980 1.2t steel sheet, in wt% Material

Mn

Si

C

P

S

DP980 1.2t

2.4800

0.1560

0.1120

0.0120

0.0044

a

Analyzed by EDS(Energy Dispersive X-ray Spectrometer)

Table 5 Mechanical properties of the DP980 1.2t at the three loading directions Loading

YS(a) [MPa]

UTS [MPa]

Elu [%]

Elf [%]

rθ (b)

rθ (c)

0

704

1014

4.44

11.2

0.835

0.840

45

690

986

4.71

12.1

1.109

1.126

90

697

1030

4.46

11.5

0.936

0.963

direction [°]

(a)

0.2% offset

(b)

From the conventional r-value calculation (in an averge sense)

(c)

From the instantaneous r-value after saturated

Table 6 Parameters of the mixed Swift–Voce model of the DP980 1.2t steel sheet at the three loading directions Loading direction [°]

K

ε0

n

A

B

C

α

0

4193.6

0.00060

0.09228

304.57

164.07

213.62

0.2226

45

3732.1

0.00074

0.09246

295.90

178.69

213.03

0.2424

90

4251.5

0.00072

0.09717

287.97

213.80

254.68

0.2277

Table 7 Equivalent plastic strains to fracture evaluated at the three loading directions Loading direction [°]

0

45

90

Uniaxial tension

Pure shear

Plane strain tension

Test#1

0.509

0.934

0.279

Test#2

0.545

0.962

0.313

Test#3

0.527

0.941

0.294

Test#1

0.640

0.783

0.255

Test#2

0.691

0.822

0.269

Test#3

0.672

0.806

0.262

Test#1

0.565

0.847

0.215

Test#2

0.603

0.887

0.246

Test#3

0.592

0.877

0.234

107

Table 8 Typical anisotropy coefficients for various materials [76] rn

∆r

1.30 to 2.00

up to 0.70

Stainless steel

0.70 to 1.10

−0.25 to 0.20

TRIP steels

0.90

about −0.03

Aluminum alloys

0.60 to 0.80

−0.60 to −0.15

Copper

0.60 to 0.80

-

Brass

0.60 to 1.00

-

Zinc alloys

0.20 to 0.60

-

Titanium alloys

2.00 to 8.00

up to 4.00

Material Deep drawing steels (DC01–DC07; cold rolled)

Table 9 Parameters of the asymptotic regression model at the three loading directions Loading direction [°]







0

0.8397

0.2310

203.37

45

1.1255

0.4239

174.51

90

0.9631

0.2266

145.06

Table 10 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the r-values obtained in an average sense F

G

H

L

M

N

0.4862

0.5450

0.4550

1.5000

1.5000

1.6591

Table 11 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the instantaneous r-values after saturated F

G

H

L

M

N

0.4741

0.5435

0.4565

1.5

1.5

1.6545

108

¯ ξ, χ and LP (σ1 ≥ σ2 ≥ σ3 ) Table 12 Typical stress condition with corresponding values of θ, θ, Typical stress condition with a hydrostatic pressure    (σ1 − σ2 ) (1, 0, 0) + σ2 (1, 1, 1) Uniaxial tension Hydrostatic pressure (σ1 , σ2 , σ2 ) =   (σ2 − σ1 ) (0, 1, 1) + σ1 (1, 1, 1) Equi-biaxial compression

θ

θ¯

ξ

χ

LP −1

0

1

1

0

π 6

0

0

0.5

0

π 3

−1

−1

1

1

Hydrostatic pressure

 σ1 + σ3 σ1 − σ3   (1, 0, −1) + (1, 1, 1)   2 2   In-plane shear Hydrostatic pressure          1  σ1 + σ3 (σ1 − σ3 ) 1, , 0 + σ3 (1, 1, 1) σ1 , , σ3 = 2 Hydrostatic pressure  2  Plane strain tenson       1   (σ3 − σ1 ) 0, , 1 + σ1 (1, 1, 1)    2 Hydrostatic pressure  Plane strain compression

(σ1 , σ1 , σ3 ) =

   (σ1 − σ3 ) (1, 1, 0) +

σ3 (1, 1, 1) Hydrostatic pressure

  (σ3 − σ1 ) (0, 0, 1) +

σ1 (1, 1, 1)

Equi-biaxial tension

Uniaxial compression

Hydrostatic pressure

Table 13 Combination of the stress triaxiality and the Lode parameter denoting typical stress states Loading condition

ηv

LP

Equi-biaxial compression

−0.6667

−1

Plane strain compression

−0.5774

0

Uniaxial compression

−0.3333

1

Pure shear

0

Uniaxial tension

0.3333

−1

Plane strain tension

0.5774

0

Equi-biaxial tension

0.6667

1

109

0

Table 14 Combination of the anisotropic stress triaxiality and the Lode parameter for the DP980 1.2t denoting the typical stress states where the parameters of Hill’s 48 criterion are obtained by the r-value evaluated in an average sense ηH

Loading condition

LP

θS = 0°

θS = 45°

θS = 90°

Equi-biaxial compression

−0.6667

−0.6516

−0.6608

−1

Plane strain compression

−0.5793

−0.5869

−0.5880

0

Uniaxial compression

−0.3333

−0.3460

−0.3460

1

Pure shear

0

0

0

0

Uniaxial tension

0.3333

0.3258

0.3304

−1

Plane strain tension

0.5793

0.5869

0.5880

0

Equi-biaxial tension

0.6667

0.6920

0.6920

1

Table 15 Model parameters of the anisotropic ductile fracture criterion C1

C2

C3

C4

7.498

2.035

1.350

1/3

Table 16 Model parameters of the anisotropic ductile fracture criterion constructed in a purely empirical way i

Pi

Qi

Ri

1

6.063

7.679

8.437

2

2.216

2.149

2.790

3

0.944

0.944

0.944

Table 17 Anisotropic characteristic parameters of the Hill’s 48 criterion determined by the normalized yield stresses where the equi-biaxial yield stress is set as the yield stress for the rolling direction F

G

H

L

M

N

0.510

0.490

0.510

1.500

1.500

1.582

Table 18 Conditions for the hydraulic bulge test of the DP980 1.2t Specimen size [mm]

Punch speed [mm/s]

Clamping force [kN]

Strain measurement

200 x 200

10

800

3D-DIC (ARAMIS)

110

Table 19 Parameters of the mixed Swift–Voce model of the DP980 1.2t steel sheet in the large deformation range Loading direction [°]

K

ε0

n

A

B

C

α

0

1925.7

0

0.3107

855.06

290.22

132.43

0.1727

45

1795.3

0

0.3764

837.88

310.27

155.24

0.1743

90

1891.4

0

0.4973

829.56

402.11

149.40

0.1535

Table 20 Mechanical properties of the DP980 1.2t steel sheet evaluated from the hydraulic bulge test YS(a) [MPa]

rb (b)

706

1.041

(a)

0.2% offset

(b)

rb = εyy /εxx

Equivalent plastic strain at the onset of fracture Test#1

Test#2

Test#3

0.616

0.723

0.669

Table 21 Anisotropic characteristic parameters of the Hill’s 48 criterion calibrated by the saturated values of the normalized yield stress ratios F

G

H

L

M

N

0.476

0.524

0.507

1.5000

1.5000

1.666

Table 22 Model parameters of the anisotropic ductile fracture criterion constructed under the proportional loading condition i

Pi

Qi

Ri

CO

1

2.901

4.616

4.857

2

2.900

2.542

3.098

3

0.669

0.669

0.669

4

0.553

0.262

0.425

1/3

Table 23 Equivalent plastic strains to fracture evaluated at the additional loading directions Loading direction [°]

22.5

67.5

Uniaxial tension

Pure shear

Plane strain tension

Test#1

0.601

0.863

0.272

Test#2

0.570

0.897

0.279

Test#3

0.620

0.889

0.283

Test#1

0.602

0.847

0.249

Test#2

0.637

0.826

0.244

Test#3

0.623

0.863

0.231

111

Table 24 Numerical algorithm for the tangent cutting plane integration method 1) Initialize k = 0 (0)

p(0)

p(0)

(0)

(0)

˜ (n+1) = σ ˜ (n) + C∆ε˜ ∆λ(n+1) = 0, ε¯(n+1) = ε¯(n) , σ 2) Check the yield criterion     (k) p(k) (k) ˜ If F(n+1) = fy σ(n+1) − h ε¯(n+1) ≥ 0 2-1) Compute the increment of the plastic multiplier factor ∂fy (k) ˜ (n+1) = m ˜ σ˜ (k) ∂σ (n+1) ∂fp (k) ˜ (n+1) = n ˜ σ˜ (k) ∂σ (n+1) dh (k) H(n+1) = d¯ εp ε¯p(k) (n+1) fp (k) R(n+1) = fy σ˜ (k) (n+1)     (k) p(k) ˜ (n+1) − h ε¯(n+1) fy σ (k) δλ(n+1) = (k) (k) (k) (k) m ˜ (n+1) : C : n ˜ (n+1) + H(n+1) R(n+1) 2-2) Update the state variables (k+1)

(k)

(k)

∆λ(n+1) = ∆λ(n+1) + δλ(n+1) p(k+1)

p(k)

(k+1)

(k+1)

(k)

(k+1)

(k)

ε¯(n+1) = ε¯(n+1) + ∆λ(n+1) R(n+1) (k)

˜ (n+1) = σ ˜ (n+1) − ∆λ(n+1) C : n ˜ (n+1) σ k = k + 1 and go back to 2) Else (k)

p(k)

˜ (n+1) = σ ˜ (n+1) , ε¯(n+1) = ε¯(n+1) σ End if 3) Go to next step

112

Table 25 Conditions for the square cup drawing test of the DP980 1.2t steel sheet Specimen size [mm]

Punch speed [mm/s]

Clamping force [kN]

299 x 299

10

1,400

Table 26 Specifications of 3,000 kN hydraulic servo press Maximum capacity [kN] and Stroke [mm]

Size [mm3 ] 3,500 x 4,400 x 4,600

Inner cylinder

Outer cylinder

Cushion cylinder

3,000, 750

1,000, 750

2,000, 350

Table 27 Total number of elements and nodes generated for each FE model Quarter FE model

Mesh size [mm]

Half FE model

Node

Element

Node

Element

3

18,207

15,000

36,057

30,000

2

40,432

33,750

80,332

67,500

1

159,607

135,000

317,100

269,100

0.5

630,000

536,406

1,257,900

1,072,812

Table 28 Parameters of the exponential fitting function for each simulation group Group

P0

P1

P2

Quarter FE model

18.763

0.160

0.000178

Half FE model

20.333

1.468

0.000248

Table 29 Anisotropic characteristic parameters of the Yld91 criterion determined by the normalized yield stress ratios saturated c1

c2

c3

c4

c5

c6

1.010

0.990

1.004

1.000

1.000

1.023

113

(a)

(b) Figure 1 Schematic forming limits classified by: (a) Necking, wrinkling and shear fracture (after Martins et al. [127]); (b) Deformation modes of the circular grid

114

Figure 2 Fracture of Dual-Phase (DP) steel sheets in various sheet metal forming processes (after research works in ductile fracture of advanced high-strength steel sheets [1, 112, 124, 155, 179, 190])

115

ˇ Figure 3 Nucleation, growth, and coalescence of voids during the tension test (after Sebek [154])

116

(a)

(b)

(c) Figure 4 Drawings of test specimens: (a) Pure shear specimen; (b) Dog-bone specimen; (c) Plane strain grooved specimen [mm]

117

Figure 5 Uniaxial tensile testing machine (INSTRON5583)

1200

RD

1100

True stress [MPa]

True stress [MPa]

1200

1000 900 800 Experiment

700

Voce fit Swift fit Mixed Swift

0

0.01

-Voce fitting curve 0.02

DD

1100 1000 900 800 Experiment

700

Voce fit Swift fit Mixed Swift

0.2226

0.03

0.04

0

True plastic strain

0.01

1100

True stress [MPa]

True stress [MPa]

1200

1000 900 800 Experiment Voce fit Swift fit Mixed Swift

0

0.01

-Voce fitting curve 0.02

0.03

0.04

(b)

TD

700

0.02

0.2424

True plastic strain

(a) 1200

-Voce fitting curve

0.03

DP980 1.2t

1100 1000 900 800

Experiment RD

700

DD TD Mixed Swift

0.2277

0.04

0

0.01

-Voce fitting curve 0.02

0.03

0.04

True plastic strain

True plastic strain

(c)

(d)

Figure 6 True stress–true strain curves of the DP980 1.2t steel sheet evaluated at the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) Flow stress curves plotted all together

118

0

Test #1 Test #2 Test #3

-1

0

Test #1 Test #2 Test #3

-1

-2 0.0

1

Pure shear (DD)

0.4

0.6

0.8

1.0

1.2

0.0

Equivalent plastic strain

Pure shear (TD)

0

Test #1 Test #2 Test #3

-1

-2 0.2

Strain path,

1

Pure shear (RD)

Strain path,

Strain path,

1

-2 0.2

0.4

0.6

0.8

1.0

1.2

0.0

Equivalent plastic strain

0.2

0.4

0.6

0.8

1.0

1.2

Equivalent plastic strain

(a)

0

Test #1 Test #2 Test #3

-1

0

Test #1 Test #2 Test #3

-1

-2 0.0

1

Uniaxial tension (DD)

0.4

0.6

0.8

1.0

0.0

Equivalent plastic strain

Uniaxial tension (TD)

0

Test #1 Test #2 Test #3

-1

-2 0.2

Strain path,

1

Uniaxial tension (RD)

Strain path,

Strain path,

1

-2 0.2

0.4

0.6

0.8

1.0

0.0

Equivalent plastic strain

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(b)

0

Test #1 Test #2 Test #3

-1

0

Test #1 Test #2 Test #3

-1

-2 0.0

1

Plane strain (DD)

0.4

0.6

0.8

Equivalent plastic strain

1.0

0.0

Plane strain (TD)

0

Test #1 Test #2 Test #3

-1

-2 0.2

Strain path,

1

Plane strain (RD)

Strain path,

Strain path,

1

-2 0.2

0.4

0.6

0.8

Equivalent plastic strain

1.0

0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(c) Figure 7 Strain paths evaluated at the three loading directions (0°, 45°, 90°) under typical loading conditions: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case

119

DP980 1.2t (Pure shear)

0.8

0.9

Fracture strain

Fracture strain

1.0

0.8

0.7

DP980 1.2t (Uniaxial tension)

0.7

0.6

0.5

Test

Test

Average value

Average value

0.0

0.0 0 (RD)

/4 (DD)

Loading direction,

0 (RD)

/2 (TD)

/4 (DD)

s

(a) 0.4

/2 (TD)

Loading direction,

s

(b)

DP980 1.2t (Plane strain tension)

1.0

Average value

Fracture strain

Fracture strain

0.8 0.3

0.2

0.6

0.4

0.2

Test

Loading condition

Average value

Pure shear

Uniaxial tension

Plane strain tension

0.0

0.0 0 (RD)

/4 (DD)

Loading direction,

/2 (TD)

0 (RD)

/4 (DD)

Loading direction,

s

(c)

/2 (TD)

s

(d)

Figure 8 Equivalent plastic strains at the onset of fracture evaluated at the three loading directions (0°, 45°, 90°) under typical loading conditions: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case; (d) Comparison among the ones averaged

120

Strain paths for each loading direction

1 RD

Strain path,

DD 0

Plane strain

TD Uniaxial tension

Pure shear

-1

-2 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

Figure 9 Representative strain paths evaluated at the three loading directions (0°, 45°, 90°) with typical loading conditions

121

0.01

RD

Transverse strain

Test #1 Test #2

0.00

Test #3 Test #4

-0.01

-0.02

r0(avg.)=0.835 -0.03 0.00

0.02

0.04

Longitudinal strain

(a) 0.01

DD

Transverse strain

Test #1 Test #2

0.00

Test #3 Test #4

-0.01

-0.02

r45(avg.)=1.109 -0.03 0.00

0.02

0.04

Longitudinal strain

(b) 0.01

TD

Transverse strain

Test #1 Test #2

0.00

Test #3 Test #4

-0.01

-0.02

r90(avg.)=0.936 -0.03 0.00

0.02

0.04

Longitudinal strain

(c) Figure 10 Transverse strain versus longitudinal strain distribution of specimens at the three loading directions: (a) 0°; (b) 45°; (c) 90°

122

RD

0.2 R

2

RD Residual

0.6642

Exp. data

Fitted data

1.2

Residuals

Lankford coefficient

1.5

0.9

0.6

0.3

0.1

0.0

Experiment data From the curve fit From the conventional

r-value

calculation

-0.1

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

(a) DD

0.2 R

2

DD Residual

0.9260

Exp. data

Fitted data

1.2

Residuals

Lankford coefficient

1.5

0.9

0.6

0.3

0.1

0.0

Experiment data From the curve fit From the conventional

r-value

calculation

-0.1

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

(b) TD

0.2 R

2

TD Residual

0.6554

Exp. data

Fitted data

1.2

Residuals

Lankford coefficient

1.5

0.9

0.6

0.3

0.1

0.0

Experiment data From the curve fit From the conventional

r-value

calculation

-0.1

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

(c) Figure 11 Instantaneous r-values fitted by the asymptotic regression model and the corresponding residuals: (a) 0°; (b) 45°; (c) 90°

123

RD

0.0004

Experiment data

-0.005

RD Residual

Exp. data

Fitted data

From the curve fit

Residual

Transverse plastic strain

0.000

-0.010

-0.015

0.0002

0.0000

-0.0002

-0.020

-0.0004

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

0.000

DD

0.0004

Experiment data

-0.005

DD Residual

Exp. data

Fitted data

From the curve fit

Residual

Transverse plastic strain

(a)

-0.010

-0.015

0.0002

0.0000

-0.0002

-0.020

-0.0004

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

0.000

TD

0.0004

Experiment data

-0.005

TD Residual

Exp. data

Fitted data

From the curve fit

Residual

Transverse plastic strain

(b)

-0.010

-0.015

0.0002

0.0000

-0.0002

-0.020

-0.0004

0.00

0.01

0.02

0.03

0.04

0.05

0.00

Longitudinal plastic strain

0.01

0.02

0.03

0.04

0.05

Longitudinal plastic strain

(c) Figure 12 Transverse plastic strain versus longitudinal plastic strain distribution predicted by the intantaneous r-values fitted by the asymptotic regression model and the corresponding residuals: (a) 0°; (b) 45°; (c) 90°

124

Equivalent plastic strain 0.01

0.02

0.03

0.04 0.00

1.0

-0.03

Normal anisotropy

0.9

-0.06

Normal anisotropy evaluated by a linear fit Planar anisotropy Planar anisotropy evaluated by a linear fit

0.8

-0.09

0.7

-0.12

0

Planar anisotropy

Normal anisotropy

0 1.1

-0.15 0

5

10

15

20

25

30

35

40

Plastic work per unit volume [MPa]

Figure 13 Instantaneous normal and plannar anisotropic quantities according to the plastic work and the equivalent plastic strain DP980 1.2t

Normalized yield stress (

yy

/

0

)

1.5

1.0

0.5 Hill's 48 (r-values averaged ) Hill's 48 (r-values saturated ) von Mises

0.0 0.5

1.0

Normalized yield stress (

1.5

xx

/

)

0

Figure 14 Comparison of the normalized yield loci constructed with Hill’s 48 and von Mises criteria 125

(a)

(b) Figure 15 Representation of the principal stress state in Haigh–Westergaard space (σ1 , σ2 , σ3 ) and on the deviatoric plane (s1 , s2 , s3 )

126

(a)

(b) Figure 16 Dependence of: (a) The stress triaxiality on various indicators concerning the deviatoric stress state; (b) Various state variables on the normalized third invariant of the deviatoric stress tensor

127

Figure 17 Rotation axis and angles for the coordinate transformation

128

(a)

(b)

129

(c)

(d) Figure 18 Stress states in the space of the newly defined stress triaxiality ηH and the Lode parameter LP according to the three loading directions of 0°, 45°, and 90°: (a) The change of new combination of the ηH and LP for DP980 1.2t steel sheet; (b) The change of new combination of the ηH and LP for SPCC 1.0t steel sheet; (c) The ηH corresponding to the typical stress paths for SPCC 1.0t steel sheet; (d) The LP corresponding to the typical stress paths 130

(a)

(b) Figure 19 Cut-off region: (a) The initial one; (b) The modified one 131

(a)

(b)

(c) Figure 20 Strain-based 3D fracture envelopes according to the three loading directions: (a) 0°; (b) 45°; (c) 90° 132

Uniaxial

Pure shear

compression

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

2D fracture locus

1.5 Loading path Pure shear Uniaxial tension Plane strain

1.0

0.5

0.0 -0.25

0.00

0.25

Stress triaxiality,

0.50

0.75

Equivalent plastic strain to fracture

Equivalent plastic strain to fracture

Fracture locus of DP980 1.2t, RD 2.0

Fracture locus of DP980 1.2t, DD 2.0

Pure shear

Uniaxial compression

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

2D fracture locus

1.5 Loading path Pure shear Uniaxial tension Plane strain

1.0

0.5

0.0 -0.25

0.00

0.25

Stress triaxiality,

H

(a)

0.50

0.75

H

(b)

2.0

Uniaxial

Pure shear

compression

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

2D fracture locus

1.5 Loading path Pure shear Uniaxial tension Plane strain

1.0

0.5

0.0 -0.25

0.00

0.25

Stress triaxiality,

0.50

Equivalent plastic strain to fracture

Equivalent plastic strain to fracture

Fracture loci of DP980 1.2t Fracture locus of DP980 1.2t, TD

2.0

Pure shear

Uniaxial compression

Plane strain

Equi-biaxial

tension

tension

tension

Loading direction o

0

1.5

o

90 (TD)

1.0 Fracture locus 0

o

45 90

(RD) o

o

(DD) (TD)

0.5

0.0 0.00

0.25

Stress triaxiality,

H

(c)

(RD) o

45 (DD)

-0.25

0.75

Uniaxial

0.50

0.75

H

(d)

Figure 21 Fracture loci according to the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) Fracture loci plotted all together

133

Damage evolution (RD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(a)

Damage evolution (DD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(b)

Damage evolution (TD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(c) Figure 22 Damage evolutions evaluated according to the loading paths and the three loading directions: (a) RD; (b) DD; (c) TD 134

FFLD of DP980 1.2t, RD

tension

tension

FFLD of DP980 1.2t, DD

2.0

Uniaxial

Plane strain

tension

tension

Major strain,

Major strain,

2.0 1

Plane strain

1

Uniaxial

1.5

Pure shear

1.5

Equi-biaxial

1.0

tension

Pure shear

Equi-biaxial

1.0

tension

Uniaxial Uniaxial

compression

FFLD

0.5

compression

FFLD

0.5

Strain path

Strain path From tests From tests

Pure shear

Pure shear

Uniaxial tension

Uniaxial tension

Plane strain

Plane strain

-1.0

-0.5

0.0 Minor strain,

0.5

1.0

-1.0

-0.5

0.0 Minor strain,

2

(a)

tension

tension

FFLDs of DP980 1.2t

2.0

Uniaxial

Plane strain

tension

tension

2.0

Loading direction o

Major strain,

Major strain,

1

1

Plane strain

1.0

2

(b)

FFLD of DP980 1.2t, TD Uniaxial

0.5

1.5

0

(RD) o

45 (DD) 90 1.5

o

(TD)

Principal stress direction [ 0

o

s

]

(RD) o

45 (DD) o

90 (TD) Pure shear

Equi-biaxial

1.0

Strain paths

Pure shear

1.0

tension

Equi-biaxial tension

Uniaxial

Uniaxial

compression

compression

FFLD

0.5

0.5

Strain path

From tests Pure shear Uniaxial tension Plane strain

-1.0

-0.5

0.0 Minor strain,

0.5

1.0

-1.0

-0.5

0.0

0.5

Minor strain,

2

(c)

1.0

2

(d)

Figure 23 Fracture forming limit diagram (FFLD) according to the three loading directions: (a) 0°; (b) 45°; (c) 90°; (d) FFLDs plotted all together

135

Equivalent plastic strain to fracture

Fracture loci of DP980 1.2t 2.0

Pure shear

Uniaxial compression

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

Loading direction 0

1.5

o

(RD) o

45 (DD) o

90 (TD)

1.0 Fracture locus 0

o

15 30

0.5

(RD) o

o

A

o

45 (DD) 60 75

A

o

o

o

90 (TD)

0.0 -0.25

0.00

0.25

0.50

Stress triaxiality,

0.75

H

(a)

FFLDs of DP980 1.2t Plane strain

tension

tension

2.0

Principal stress direction [ 0

Major strain,

1

Uniaxial

o

]

s

(RD)

15 30

o

o

o

45 (DD) 60

1.5

75

o

o

o

90 (TD)

Loading direction o

0

(RD)

Strain paths

o

45 (DD) o

Pure shear

90 (TD)

1.0

Equi-biaxial tension

A

Uniaxial compression

0.5 B

B

A

-1.0

-0.5

0.0 Minor strain,

0.5

1.0

2

(b) Figure 24 Fracture loci and its corresponding FFLDs at different directions of the maximum principal stress: (a) Fracture loci plotted according to the maximum principal stress direction; (b) FFLDs plotted according to the maximum principal stress direction

136

Major fracture strain in plane strain 0.28

FFLD0

0.26

0.24 FFLD

FFLD

0

0

predicted

from the fracture criterion

RD

0.22

DD TD

0.00 0 (RD)

/4 (DD)

Max. principal stress direction,

/2 (TD)

s

Figure 25 Major fracture strain in plane strain according to the maximum principal stress direction

Figure 26 Scaling factor in the space of the maximum principal stress direction and the strain path

137

FFLDs of DP980 1.2t scaled Plane strain

tension

tension

2.0

Major strain,

1

Uniaxial

1.5 FFLD (Reference ) o

0

(RD)

Strain path

Pure shear

FFLD (scaled)

1.0

45 90

o

o

Equi-biaxial tension

(DD) (TD)

Uniaxial compression

0.5

-1.0

-0.5

0.0

0.5

Minor strain,

1.0

2

(a)

Comparison of the FFLDs with the scaled ones Plane strain

tension

tension

2.0

Principal stress direction [ 0

Major strain,

1

Uniaxial

o

]

s

(RD)

15 30

o

o

o

45 (DD) 60

1.5

75

o

o

o

90 (TD) Scaled from the FFLD for the RD o

45 (DD) o

90 (TD)

Strain paths

Loading direction Pure shear

o

0

1.0

(RD) o

45 (DD)

Equi-biaxial

o

90 (TD)

tension

A

Uniaxial compression

0.5 B

B

A

-1.0

-0.5

0.0 Minor strain,

0.5

1.0

2

(b) Figure 27 Fracture forming limit diagrams at different directions of the maximum principal stress: (a) FFLDs scaled at different directions of the maximum principal stress; (b) Comparison between FFLDs predicted from the anisotropic ductile fracture criterion and the FFLDs scaled

138

(a)

(b)

(c) Figure 28 Stress-based 3D fracture envelopes according to the three loading directions: (a) 0°; (b) 45°; (c) 90°

139

1

Stress-based FFLDs of DP980 1.2t

Major true stress,

2000

1500

Uniaxial tension

Plane strain Equi-biaxial tension

1000

Pure shear

FFLD (Hill's 48 anisotropic ) 0

o

45 90

500

(RD) o

o

(DD) (TD)

Loading path

Yield locus von Mises Hill's 48

-1000

-500

0

500

Minor true strss,

1000

1500

2

(a)

Stress-based fracture loci of DP980 1.2t

Equivalent stress to fracture [MPa]

2400 Loading direction 0

o

(RD) o

45 (DD)

1800

o

90 (TD)

Fracture locus 0

1200

o

(RD) o

45 (DD) o

equi-biaxial

90 (TD)

compression

Yield

600

Hill's 48 Plane strain

Uniaxial compression

Uniaxial tension

compression

Plane strain tension

Pure shear

equi-biaxial tension

0 -0.75

-0.50

-0.25

0.00

Stress triaxiality,

0.25

0.50

0.75

H

(b) Figure 29 Stress-based fracture forming limit criteria according to the three loading directions (0°, 45°, 90°): (a) The stress-based FFLDs; (b) The corresponding fracture loci

140

1

Stress-based FFLDs of DP980 1.2t Major true stress,

2000

1500

Uniaxial tension

Plane strain Equi-biaxial tension

1000

Pure shear

Fracture locus 0

o

(RD)

15 500

Loading path

30

o

o

o

45 (DD) 60

Yield locus

-1000

-500

0

o

von Mises

75

Hill's 48

90 (TD)

500

Minor true strss,

o

o

1000

1500

2

(a)

Stress-based fracture loci of DP980 1.2t

Equivalent stress to fracture [MPa]

2100 Compression region

Shear

dominated region

FLD region

Loading direction 0

o

(RD) o

45 (DD)

1800

o

90 (TD)

Fracture locus 0

1500

o

15 equi-biaxial

30

compression

(RD) o

o

o

45 (DD) 60

1200

Plane strain

Uniaxial compression

Uniaxial tension

75

Plane strain tension

compression

Pure shear

o

o

o

90 (TD)

equi-biaxial tension

0 -0.75

-0.50

-0.25

0.00

Stress triaxiality,

0.25

0.50

0.75

H

(b) Figure 30 Stress-based fracture forming limit criteria according to various directions of the maximum principal stress: (a) The stress-based FFLDs; (b) The corresponding fracture loci

141

(a)

(b) Figure 31 Stress-based fracture forming limit diagram: (a) Not considering the one point convergence at the equi-biaxial fracture stress; (b) Considering the one point convergence at the equi-biaxial fracture stress

142

PEPS FFLDs of DP980 1.2t 2.0

Uniaxial

Plane strain tension

tension

Predicted from the fracture model 0

o

(RD) o

45 (DD) o

90 (TD)

1.5

Loading direction 0

o

(RD) o

45 (DD) o

90 (TD) Equi-biaxial Pure shear

tension

1.0

Strain paths

Uniaxial compression

-1.0

0.5

-0.5

0.0

0.5

1.0

Figure 32 Polar Effective Plastic Strain (PEPS) Fracture Forming Limit Diagrams (FFLDs) according to the three loading directions (0°, 45°, 90°)

143

PEPS FFLDs of DP980 1.2t 2.0

Uniaxial

Plane strain tension

tension

Predicted from the fracture model 0

o

(RD) o

45 (DD) o

90 (TD)

1.5

Scaled from the FFLD for RD o

45 (DD) o

90 (TD) Loading direction 0

o

(RD)

Equi-biaxial

o

Pure shear

45 (DD)

1.0

tension

o

90 (TD)

Strain paths

Uniaxial Uniaxial compression compression

0.5

A A

-1.0

-0.5

0.0

0.5

1.0

Figure 33 Comparison of the Polar Effective Plastic Strain (PEPS) Fracture Forming Limit Diagrams (FFLDs) predicted from the fracture criterion with the FFLDs scaled from the strainbased FFLD for RD

144

Damage evolution (RD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(a)

Damage evolution (DD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(b)

Damage evolution (TD) 1.2

Damage

1.0

0.8

0.6

0.4 Pure shear 0.2

Uniaxial tension Plane strain

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Equivalent plastic strain

(c) Figure 34 Damage evolutions according to the loading paths and three loading directions when the empirical way is employed to the newly developed fracture criterion: (a) 0°; (b) 45°; (c) 90° 145

Fracture Forming Limit Diagram, RD Uniaxial

Plane strain

tension

tension

Fracture Forming Limit Diagram, DD

1.5

Uniaxial

Plane strain

tension

tension 1

From tests

Major strain,

Major strain,

Pure shear

Pure shear Uniaxial tension Plane strain

1.0

Equi-biaxial

Pure shear

tension

Uniaxial

Uniaxial

compression

compression

0.5

-1.0

-0.5

1.5 FFLD

1

FFLD

From tests Pure shear Uniaxial tension

0.0

0.5

1.0

-1.0

-0.5

0.0

tension

FFLDs of DP980 1.2t improved

1.5

Uniaxial

Plane strain

tension

tension

1.5

Principal stress direction [

1

1

FFLD

0 Major strain,

Major strain,

Pure shear

1.0

2

(b)

Fracture Forming Limit Diagram, TD Plane strain

0.5

Minor strain,

2

(a)

tension

tension

0.5

Minor strain,

Uniaxial

Equi-biaxial

Plane strain

1.0

From tests Pure shear Uniaxial tension Plane strain

Equi-biaxial tension

1.0

Pure shear

o

s

]

(RD) o

45 (DD) o

90 (TD)

Strain paths

1.0 Equi-biaxial tension

Uniaxial

Uniaxial

compression

compression

0.5

0.5

Loading direction o

0

(RD) o

45 (DD) 90

-1.0

-0.5

0.0 Minor strain,

0.5

1.0

-1.0

-0.5

0.0

(c)

(TD)

0.5

Minor strain,

2

o

2

(d)

Figure 35 Fracture forming limit diagram (FFLD) according to the three loading directions when the empirical way is employed to the newly developed fracture criterion: (a) 0°; (b) 45°; (c) 90°; (d) FFLDs plotted all together

146

1.0

Figure 36 Schematic description of the radius of the curvature

Figure 37 Hydraulic bulge test system in POSCO [7]

147

Flow stress [MPa]

1500

DP980 1.2t

1000

Hydraulic bulge test

500

Test #1 Test #2 Test #3 Test #4 Test #5

0

0.2

0.4

0.6

0.8

Equivalent plastic strain

Figure 38 Flow stress curves obtained from the hydraulic bulge tests

148

1.0

(a)

(b) Figure 39 Hydraulic bulge test result: (a) Fractured specimen; (b) Equivalent plastic strain distribution before the onset of fracture

149

Figure 40 Schematic description of the plastic work equivalence principle

150

DP980 1.2t

3

Plastic work [mJ/mm ]

40

30

20 Uniaxial tension RD DD

10

TD Equi-biaxial tension

0

0.01

0.02

0.03

0.04

0.05

Equivalent plastic strain

Normalized yield stress,

u/ b

(a)

1.10

DP980 1.2t

1.05 1.00 0.95 0.90

After saturation 0

0.85

b

0

5

0.985,

45

RD 0.961,

b

10

90

DD 1.009

TD

b

15

20

25

30

35

3

Plastic work [mJ/mm ]

(b) Figure 41 Plastic works per volume of the DP980 1.2t steel sheet and the directional yield stresses normalized by the equi-biaxial yield stress: (a) Plastic works per volume; (b) Normalized yield stress ratios

151

Flow stress [MPa]

1500

DP980 1.2t

1000 Hardening curve obtained from the combination of the tensile test and the hydraulic bulge test RD DD

500

TD Mixed Swift

Voce fitting curve

RD DD TD

0

0.1

0.2

0.3

Equivalent plastic strain

Figure 42 Flow stress curves in the large deformation range

152

0.4

(a)

(b)

(c) Figure 43 Strain-based 3D fracture envelopes constructed under the proportional loading condition: (a) 0°; (b) 45°; (c) 90° 153

Equivalent plastic strain to fracture

Fracture loci of DP980 1.2t Uniaxial compression

Pure shear

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

Loading direction o

0

b f

(RD) o

45 (DD) o

90 (TD)

Fixed

Experimental data 0

o

(RD)

45 90

o

o

(DD) (TD)

Equi-biaxial tension

Fracture locus 0

o

45 90

Stress triaxiality,

(RD) o

o

(DD) (TD)

H

(a)

(b) Figure 44 Fracture loci constructed under the proportional loading condition: (a) Predicted from the fracture criterion considering the model parameters dependent on the loading direction; (b) Predicted from the fracture criterion considering the model parameters dependent on the loading direction in combination with the weight function

154

Pure shear

Uniaxial tension

Uniaxial

Major strain,

1

FFLDs of DP980 1.2t 1.0

Plane strain

Equi-biaxial

tension

tension

Loading direction

0.8

0

o

(RD)

45 90

0.6

o

o

Strain paths

(DD) (TD)

Experimental data

compression

0

o

(RD)

45

0.4

90

o

(DD)

o

(TD) Principal stress direction 0

0.2

o

45 90

-1.0

-0.5

0.0

s

]

(RD) o

o

(DD) (TD)

0.5

Minor strain,

[

1.0

2

(a)

Pure shear

Uniaxial tension

Uniaxial

Major strain,

1

FFLDs of DP980 1.2t improved 1.0

Plane strain

Equi-biaxial

tension

tension

Loading direction

0.8

0

o

(RD)

45 90

0.6

o

o

Strain paths

(DD) (TD)

Experimental data

compression

0

o

45

0.4

90

(RD) o

o

(DD) (TD) Principal stress direction 0

0.2

o

45 90

-1.0

-0.5

0.0 Minor strain,

[

s

]

(RD) o

o

(DD) (TD)

0.5

1.0

2

(b) Figure 45 Strain-based FFLDs constructed under the proportional loading condition predicted from the fracture criterion considering the model parameters dependent on: (a) The loading direction; (b) The loading direction in combination with the weight function

155

Pure shear

Uniaxial tension

Uniaxial

Major strain,

1

FFLDs of DP980 1.2t 1.0

Plane strain

Equi-biaxial

tension

tension

Loading direction

0.8

0

o

(RD)

45 90

0.6

o

o

Strain paths

(DD) (TD)

Experimental data

compression

0

o

45

0.4

90

(RD) o

o

(DD)

FLD of DP980 1.2t (TD)

(TD) Principal stress direction

From the isotropic Lou

-Huh ductile fracture

0

0.2

45

criterion

-1.0

o

90

-0.5

0.0 Minor strain,

0.5

[

S]

(RD) o

o

(DD) (TD)

1.0

2

(a)

(b) Figure 46 Comparison between anisotropic and isotropic fracture forming limits: (a) strainbased FFLDs with the necking-based FLD; (b) fracture loci

156

(a)

(b) Figure 47 Fracture forming limit envelopes constructed under the proportional loading condition: (a) Strain-based fracture forming envelope; (b) Stress-based fracture forming envelope

157

(a)

(b) Figure 48 Three-dimensional fracture envelopes constructed under the proportional loading condition: (a) In the principal strain space; (b) In the principal stress space

158

Equivalent plastic strain to fracture

1.0

Pure shear

Additional tests

0.9

0.8 Test Average value

0.7

Predicted from the proposed fracture model

0.0 0 (RD)

/8

/4 (DD)

/8

Loading direction,

/2 (TD)

S

Equivalent plastic strain to fracture

(a) 0.8

Uniaxial tension Test Additional tests

Average value

0.7

0.6

0.5 Predicted from the proposed fracture model

0.0 0 (RD)

/8

/4 (DD)

/8

Loading direction,

/2 (TD)

S

Equivalent plastic strain to fracture

(b) 0.4

Plane strain tension Test Average value

Additional tests

0.3

0.2

Predicted from the proposed fracture model

0.0 0 (RD)

/8

/4 (DD)

Loading direction,

/8

/2 (TD)

S

(c) Figure 49 Comparison of the fracture strains obtained from the experiments with the ones predicted from the proposed fracture criterion under the proportional loading condition: (a) Pure shear; (b) Uniaxial tension; (c) Plane strain tension 159

Figure 50 Prediction of the fracture strain over a wide range of stress state under the plane stress condition

Figure 51 Normal directions of the yield surface and the plastic potential in the non-AFR approach

160

Figure 52 Relationship between the acceleration, velocity, and the displacement in the explicit time integration scheme

Figure 53 Graphical description of the tangent cutting plane algorithm

161

Figure 54 Specimen fractured during the square cup drawing test

162

(a)

(b) Figure 55 Schematic descriptions of each square cup drawing simulation: (a) Symmetric case; (b) Non-symmetric case with a bias of 5 mm

163

DP980 1.2t

Normalized yield stress (

yy

/

0

)

1.5

1.0

0.5 Hill's 48 plastic potential Hill's 48 yield function von Mises yield function

0.0 0.5

1.0

Normalized yield stress (

1.5

xx

/

)

0

Figure 56 Yield loci and plastic potential of the DP980 1.2t steel sheet

Flow stress [MPa]

1500

DP980 1.2t

1000

Stress

500

-Strain response in large plastic strain range

(A) evaluated from the tensile test up to UTS (B) extrapolated from the tensile test result using the mixed Swfit

-Voce hardening law

(C) obtained from the bulge test result based on the plastic work equivalence

0

0.1

0.2

0.3

Equivalent plastic strain

Figure 57 Flow stress curves used in the FE analysis

164

0.4

(a)

(b)

(c) Figure 58 Finite elemement modeling of each specimen: (a) Pure shear; (b) Uniaxial tension; (c) Plane strain tension

165

Load [kN]

6

Pure shear

4

Test #1 Test #2

2

Test #3 Simulation (A) Simulation (B) Fracture simulation

0

0.5

1.0

1.5

2.0

2.5

3.0

2.0

2.5

3.0

L [mm]

Displacement,

(a)

Load [kN]

10

Uniaxial tension

8

6 Test #1

4

Test #2 Test #3 Simulation (A)

2

Simulation (B) Fracture simulation

0

0.5

1.0

1.5

Displacement,

L [mm]

(b)

Load [kN]

20

Plane strain tension

15

10

Test #1 Test #2 Test #3 Simulation (A)

5

Simulation (B) Fracture simulation

0

0.1

0.2

0.3

Displacement,

0.4

L [mm]

0.5

(c) Figure 59 Comparison of the load–displacement curves for each specimen: (a) Pure shear case; (b) Uniaxial tension case; (c) Plane strain tension case

166

Figure 60 Results of the FE analysis for a quarter FE model according to various mesh size

Figure 61 Results of the FE analysis for a half FE model according to various mesh size

167

(a)

(b) Figure 62 Crack predicted in FE simulation: (a) A quarter FE model; (b) A half FE model

168

(a)

(b)

Figure 63 Localized equivalent plastic strain distribution before the fracture initiation with different mesh sizes of: (a) 0.5 mm; (b) 3 mm

169

Fracture Stroke [mm]

40

Quarter model

From the experiment Fitted curve

30

Simulation result

20

10

0

0.02

0.04

0.06

0.08

0.10

0.12

Mesh size/Punch radius

(a)

Fracture Stroke [mm]

40

Half model

From the experiment Fitted curve

30

Simulation result

20

10

0

0.02

0.04

0.06

0.08

0.10

0.12

Mesh size/Punch radius

(b) Figure 64 Punch stroke to fracture for each FE group: (a) A quarter FE model; (b) A half FE model

170

(a)

(b) Figure 65 Fracture prediction using the strain-based FFLD for an: (a) Outer layer; (b) Inner layer

171

Figure 66 Fracture prediction using the strain-based FFLD according to the loading direction

172

(a)

(b) Figure 67 Fracture prediction of an outer layer using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence

173

(a)

(b) Figure 68 Fracture prediction of an inner layer using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence

174

(a)

(b) Figure 69 Fracture prediction according to the loading direction using the stress-based FFLD: (a) Not considering the one point convergence; (b) Considering the one point convergence 175

Figure 70 Fracture prediction using the PEPS FFLD at the punch edge of the sheet

Figure 71 Fracture prediction using the PEPS FFLD at the punch radius of the sheet

176

Figure 72 Fracture prediction using the PEPS FFLD at the punch edge of the sheet where no fracture initiates

177

(a)

(b) Figure 73 Conditions on the material parameters to reduce: (a) The Yld2004-18p criterion to the Yld91 criterion; (b) The Yld91 criterion to Hosford, Hill’s 48, von Mises, and Tresca criteria 178

Figure 74 Yield loci of AA6111-T4 constructed by the generalized anisotropic stress triaxiality under the plane stress condition

Figure 75 Rotation axis and angles defined for the application of the fracture surface to sheet metal forming

179

(a)

(b)

(c) Figure 76 Distribution of the scaling values over a wide range of the stress states according to the three loading directions: (a) 0°; (b) 45°; (c) 90°

180

1.01

Scaling value

1.00

0.99

0.98

0.97 Loading direction

0.96

o

0

(RD) o

45

(DD) o

90

(TD)

0.00 -1.0

-0.5

0.0

Lode parameter,

0.5

1.0

LP

Figure 77 Scaling values according to the Lode parameter at the three loading directions

(a) Figure 78 Influence of the exponent m of the Yld91 criterion on the yield locus

181

Equivalent plastic strain to fracture

Influence of the exponent of the Yld91 criterion 1.0

Uniaxial

Pure shear

compression

Uniaxial

Plane strain

Equi-biaxial

tension

tension

tension

0.8

0.6

0.4

Exponent of the Yld91 criterion m = 8 m = 6 m = 2

0.2

Stress triaxiality m = 8 m = 6 m = 2

0.0 -0.25

0.00

0.25

0.50

Stress triaxiality,

0.75

A

(a)

Pure shear

Uniaxial tension

Uniaxial

Major strain,

1

Influence of the exponent of the Yld91 criterion 1.0

Plane strain

Equi-biaxial

tension

tension

0.8

Strain paths

0.6

compression

Exponent of the Yld91 criterion m = 8

0.4

m = 6 m = 2 Strain path

0.2

m = 8 m = 6 m = 2

-1.0

-0.5

0.0

0.5

Minor strain,

1.0

2

(b) Figure 79 Influence of the exponent m of the Yld91 criterion on the fracture forming limit: (a) Fracture loci; (b) Fracture forming limit diagrams

182