A non-iterative parameter identification procedure for the non-local Gurson–Tvergaard–Needleman model based on standardized experiments

Abstract

Damage mechanics models exhibit favorable properties such as the intrinsic influence of stress triaxiality on damage evolution and the prediction of crack initiation as well as propagation leading to structural failure. However, their application requires advanced expertise hindering the transfer of these models into industrial practice, especially since the parameter calibration is a key obstacle. In this paper, a simplified procedure is proposed for a non-local extension of the Gurson–Tvergaard–Needleman model (GTN), which is a highly accepted model for ductile failure of metals. The procedure is iteration free and requires experimental input data from only two standardized tests. The parameters are determined using look-up diagrams created on the basis of systematic simulations and made available for different material behavior covering the majority of ductile metals. Benchmark tests for three different steels are conducted to evaluate the robustness of the proposed procedure. The reliability of the GTN model is validated for all investigated materials.

1 Introduction

Ductile metals are used in different areas of engineering. This wide employment has always drawn special attention to an in-depth investigation of ductile failure to ensure safe operation and increased the component’s efficiency with respect to its weight and load bearing capability. Ductile failure is mainly prompted by the evolution of the so-called ductile damage, which is characterized by the nucleation, growth and coalescence of microvoids due to plastic deformation.

In the last decades, the modeling of ductile damage and failure received substantial interest. Numerous modeling attempts were made to tackle and overcome the different major challenges to realistically predict and simulate ductile failure at different scales. Reviews were given, e.g., by Besson (2010), Benzerga et al. (2016) or Li and Cui (2020). Among the numerous approaches, the micromechanically-based model of Gurson (1977) with its extension of Tvergaard and Needleman (1984), famously known as the GTN model, has since been characterized as the gold standard to simulate and predict ductile damage and failure. The GTN model introduces a yield function which does not only depend on the equivalent Mises stress, but additionally on the hydrostatic stress and the void volume fraction $f$, where the latter characterizes the damage state and evolves with plastic deformation.

In the original local version, the GTN model suffers (as do all other local damage models) from a spurious mesh dependency of respective finite element simulations. The damage localizes in a thin zone, determined by the element size, at the initiation of softening. In order to overcome this problem, different non-local formulations of the GTN-model and other ductile damage models have been proposed and applied successfully (Leblond et al. 1994; Jackiewicz and Kuna 2003; Jirásek and Rolshoven 2003; Reusch et al. 2003; Linse et al. 2012; Enakoutsa 2014; Soyarslan et al. 2016; Zhang et al. 2018; Nguyena et al. 2020; Bergo et al. 2021). Due to their simple and robust integrability into existing finite element codes, implicit gradient enhanced formulations (Peerlings et al. 1996) and related approaches have become the standard tool for obtaining mesh convergent results (Linse et al. 2012; Zhang et al. 2018; Nguyena et al. 2020; Samal et al. 2008; Kiefer et al. 2018; Seupel et al. 2018, 2020; Ostwald et al. 2019; Azinpour et al. 2018).

Another challenge encounters the different damage models and is manifested by the number of parameters to be identified. For instance, the extension of Tvergaard and Needleman introduced several new parameters to accurately predict the void growth kinetics and coalescence. Additionally, the strain-induced void nucleation proposed by Chu and Needleman (1980) introduces more parameters. Furthermore, the non-local extensions add extra parameters to the identification procedure, e.g. the intrinsic length. Great efforts have been made in the past on different model calibration strategies and parameter identification procedures (cf. also the review in Li and Cui (2020)). Many practical and useful conclusions resulted from the different studies, which facilitate the identification process. For example, there is a broad consensus in the literature that the yield curve of the matrix material can be determined from a smooth tensile test. Moreover, it was found that the initial and nucleable porosity, $f_0$ and $f_\textrm{n}$, respectively, should be identified from metallographic analysis, if available (e.g. by Franklin’s formula). Regarding the other parameters of the GTN model, being related to the growth and coalescence of voids, the identification strategies differ considerably. Starting with the seminal works of Tvergaard and Needleman (1984), Tvergaard (1981), it was proposed to calibrate all these parameters from numerical unit cell simulations (Kuna and Sun 1996; Steglich and Brocks 1998; Faleskog et al. 1998; Koplik and Needleman 1988). In order to improve the predictions for the void coalescence stage, it was later proposed in several studies to calibrate the critical void volume fraction $f_{\text {c}}$ from notched tensile tests, see (Jackiewicz and Kuna 2003; Zhang et al. 2021; Haušild et al. September 2002; Brocks et al. 1995; Bernauer and Brocks 2002) and references therein. In several recent studies, all GTN parameters were calibrated with the help of notched tensile tests (Chen et al. 2021; Zhu et al. 2022) or the Small-Punch Test (Abendroth and Kuna 2003; Kulawinski et al. 2021) using generic optimization algorithms, with a certain grouping of the parameters with respect to their influence on certain regimes of the load-deflection curves (Leclerc et al. 2021). A round-robin study (Bernauer and Brocks 2002), however, revealed problems on the reproducibility of the results and rather weak predictions of fracture tests if the GTN parameters are calibrated only from notched tensile tests. Some studies therefore directly employed fracture tests for calibrating the respective parameters of the GTN model (Zhang et al. 2018; Seupel et al. 2020; Gao et al. 1998; Hütter et al. 2013a). Fracture mechanics tests, as well as smooth tensile tests, have the advantage that they are standardized. From the engineering perspective, all the aforementioned strategies of parameter calibration have in common that they require iterations over relatively expensive finite element simulations.

The lack of a reliable and yet simple parameter calibration strategy has been identified as one of the obstacles for a wide-spread usage of the GTN-model (and other damage models) in engineering (Li and Cui 2020; Bernauer and Brocks 2002). The aim of the present work is to develop and validate an easy-to-use non-iterative procedure based on standardized experiments.

The outline of the present paper is as follows. After a brief introduction, a short description of the employed implicit gradient-enhanced GTN model and its numerical implementation are presented. In the following section, the simplified parameter identification procedure is introduced. Next, benchmark examples are shown for different materials and conditions. Afterwards, some recommendations and usage guidelines are given to assist in obtaining the best outcomes. Finally, the main accomplishments and conclusions of this work are summarized.

2 Non-local GTN-model (Gradient enhanced formulation)

2.1 Formulation

In the present study, a non-local modification of the continuum damage model established by Gurson (1977) and improved by Tvergaard and Needleman (1984), Tvergaard (1981) is employed, adapting an implicit gradient-enriched formulation. In general, the void volume fraction $f$ is used as a measure for ductile damage, which enters the yield function by its effective counterpart $f^*$. The yield function $\varPhi$ of the original GTN model still holds in this formulation and is given in a rate independent formulation with loading-unloading conditions as follows

$$\begin{aligned} \varPhi&= \frac{\sigma _{\text {eq}}^2}{{\bar{\sigma }}^2} {+}2\, q_1 \, f^*\cosh \left( \frac{3\, q_2 \, \sigma _{\text {h}}}{2 \, {\bar{\sigma }}} \right) {-}\left[ 1{+} \left( q_1 \, f^*\right) ^2\right] \le 0 \, ,\nonumber \\&\quad \quad \Lambda _{\text {pl}}\ge 0 \, , \quad \Lambda _{\text {pl}}\, \varPhi =0 \, , \end{aligned}$$

(1)

in terms of the hydrostatic stress $\sigma _{\text {h}}$, the Mises equivalent $\sigma _{\text {eq}}$ stress and with $q_1$, $q_2$ as the empirical model parameters of Tvergaard and Needleman. Furthermore, ${\bar{\sigma }}$ is the current yield stress of the matrix material and $\Lambda _{\text {pl}}$ corresponds to the plastic multiplier. The mentioned characteristic stress measures are defined as

$$\begin{aligned} \sigma _{\text {h}}&= \frac{1}{3} \, \text {tr}\left( \varvec{\sigma }\right) \, ,&\varvec{S}&= \varvec{\sigma }- \sigma _{\text {h}}\, {\varvec{I}} \,,&\sigma _{\text {eq}}&= \sqrt{\frac{3}{2} \, \varvec{S}: \varvec{S}} \,, \end{aligned}$$

(2)

where $\varvec{S}$ denotes the stress deviator and ${\varvec{I}}$ is the identity tensor of second order. A hypoelastic approach and associative flow rule are employed, where the stress–strain relation is based on the assumption that the deformation rate $\varvec{D}$ can be decomposed into an elastic $\varvec{D}_{\text {el}}$ and plastic $\varvec{D}_{\text {pl}}$ part, so that

$$\begin{aligned} \dot{\varvec{\sigma }}^{\text {J}}&= {\varvec{\mathcal {C}}}:\varvec{D}_{\text {el}}\, ,&\varvec{D}_{\text {el}}&= \varvec{D}- \varvec{D}_{\text {pl}}\, . \end{aligned}$$

(3)

Therein, $\dot{\varvec{\sigma }}^{\text {J}}$ and ${\varvec{\mathcal {C}}}$ correspond to the Jaumann-rate of the Cauchy-stress tensor $\varvec{\sigma }$ and to the fourth order tensor of isotropic elasticity, respectively. The associative flow rule of the plastic rate of deformation yields

$$\begin{aligned} \varvec{D}_{\text {pl}}&= \Lambda _{\text {pl}}\, \frac{\partial \varPhi }{\partial \varvec{\sigma }}={\dot{\varepsilon }}_{\text {eq}}\, \varvec{N}+{\dot{\varepsilon }}_{\text {h}}\, {\varvec{I}} ,&\varvec{N}&=\frac{3}{2 \, \sigma _{\text {eq}}} \, \varvec{S}\, , \end{aligned}$$

(4)

using the definitions of the rates of macroscopic equivalent plastic strain ${\dot{\varepsilon }}_{\text {eq}}$ and the volumetric plastic strain ${\dot{\varepsilon }}_{\text {h}}$

$$\begin{aligned} {\dot{\varepsilon }}_{\text {eq}}&= \Lambda _{\text {pl}}\, \frac{\partial \varPhi }{\partial \sigma _{\text {eq}}} \, ,&{\dot{\varepsilon }}_{\text {h}}&=\Lambda _{\text {pl}}\, \frac{1}{3} \, \frac{\partial \varPhi }{\partial \sigma _{\text {h}}} \, . \end{aligned}$$

(5)

The effective equivalent plastic strain of the matrix material ${\bar{\varepsilon }}$ is determined by the evolution equation

$$\begin{aligned} \dot{{\bar{\varepsilon }}}&= \frac{\varvec{\sigma }:\varvec{D}_{\text {pl}}}{\left[ 1-f\right] \sigma _{\text {y}}\left( {\bar{\varepsilon }}\right) } = \frac{\sigma _{\text {eq}}\, {\dot{\varepsilon }}_{\text {eq}}+3 \, \sigma _{\text {h}}\, {\dot{\varepsilon }}_{\text {h}}}{\left[ 1-f\right] \sigma _{\text {y}}\left( {\bar{\varepsilon }}\right) } \, . \end{aligned}$$

(6)

The evolution of the void volume fraction $f$ is related to two main contributions, the growth of the existing voids ${\dot{f}}_\textrm{G}$ and the nucleation of new voids ${\dot{f}}_\textrm{N}$

$$\begin{aligned} {\dot{f}}= {\dot{f}}_\textrm{G} + {\dot{f}}_\textrm{N}. \end{aligned}$$

(7)

The growth of the initial voids is related to the volumetric plastic strain evolution ${\dot{\varepsilon }}_{\text {h}}$, since an incompressible matrix behavior is assumed, which is expressed by the relation

$$\begin{aligned} {\dot{f}}_\textrm{G} = 3 \left[ 1-f\right] \, {\dot{\varepsilon }}_{\text {h}}\, . \end{aligned}$$

(8)

The void nucleation is driven by the equivalent plastic strain ${\dot{\varepsilon }}_{\text {eq}}$

$$\begin{aligned} {\dot{f}}_\textrm{N} = {\mathcal {A}}\left( \varepsilon _\textrm{eq}\right) \, {\dot{\varepsilon }}_{\text {eq}}\, , \end{aligned}$$

(9)

where the specific form of the nucleation rate ${\mathcal {A}}$ is chosen according to Chu and Needleman (1980)

$$\begin{aligned} {\mathcal {A}} \left( \varepsilon _\textrm{eq}\right) = \frac{f_\textrm{n}}{s_\textrm{n}\sqrt{2\pi }}\,\textrm{exp} \left( -\frac{1}{2}\left[ \frac{\varepsilon _\textrm{eq} -\varepsilon _\textrm{n}}{s_\textrm{n}}\right] ^2\right) \, . \end{aligned}$$

(10)

Therein, $f_\textrm{n}$ is the volume fraction of nucleable voids and $\varepsilon _\textrm{n}$ as well as $s_\textrm{n}$ correspond to the parameters of the normal distribution.

In the considered non-local modification of the GTN-model, the evolution equation of the void volume fraction is changed. A local strain-like quantity ${\dot{\varepsilon }}_{\text {l}}$ is introduced which comprises the void nucleation and void growth according to

$$\begin{aligned} {\dot{\varepsilon }}_{\text {l}}&= \frac{\dot{f}}{3 \left[ 1-f\right] } = {\dot{\varepsilon }}_{\text {h}}+ \frac{1}{3 \left[ 1-f\right] }\,{\mathcal {A}} \left( \varepsilon _\textrm{eq}\right) \, {\dot{\varepsilon }}_{\text {eq}}\, . \end{aligned}$$

(11)

The non-local counterpart $\varepsilon _{\text {nl}}$ related to $\varepsilon _{\text {l}}$ is introduced following the implicit gradient approach of Peerlings et al. (1996) by extending the mechanical boundary value problem with a further partial differential equation of Helmholtz type

$$\begin{aligned} l_{\text {nl}}^2 \, \Delta _{{x}}\varepsilon _{\text {nl}}&=\varepsilon _{\text {nl}}- \varepsilon _{\text {l}},&\forall x \in \Omega \, . \end{aligned}$$

(12)

Trivial Neumann-boundary conditions are prescribed on all surfaces

$$\begin{aligned} \text {grad}_{{x}}\varepsilon _{\text {nl}}\cdot \vec {n}&= 0 ,&\forall x \in \partial \Omega \, . \end{aligned}$$

(13)

The symbols $\Delta _{{x}}$, $\text {grad}_{{x}}$ and $\vec {n}$ are the Laplacian operator, the gradient operator and the outward unit normal vector with respect to the current configuration, respectively. Equation (12) introduces an additional parameter $l_{\text {nl}}$, which can be interpreted as an internal length that controls the width of the damage process zone.

The mechanical boundary value problem remains unchanged in the non-local formulation and is given by the balance of linear (static case) and angular momentum

$$\begin{aligned} \text {div}_{{x}}\varvec{\sigma }&={\vec {0}} \, ,&\varvec{\sigma }&={\varvec{\sigma }}^{\text {T}}&\forall x&\in \Omega \end{aligned}$$

(14)

and the respective boundary conditions.

The regularization of the boundary value problem is attained by relating the void growth ${\dot{f}}$ to the evolution of the non-local plastic strain ${\dot{\varepsilon }}_{\text {nl}}$

$$\begin{aligned} {\dot{f}}&=3 \, \left[ 1-f\right] \, {\dot{\varepsilon }}_{\text {nl}}\, . \end{aligned}$$

(15)

For vanishing void nucleation, the non-local variable $\varepsilon _{\text {nl}}$ corresponds to the non-local volumetric plastic strain as introduced by Linse et al. (2012). Moreover, a full regularization is obtained even in the case of a nucleation dominated damage process. However, in contrast to the recent GTN-approaches by Nguyena et al. (2020) and Leclerc et al. (2021), void nucleation and growth are handled by the same length scale parameter $l_{\text {nl}}$.

A modified accelerated void evolution formulation is employed here to model the coalescence mechanism (Seupel et al. 2020)

$$\begin{aligned} f^*\left( f\right)&= {\left\{ \begin{array}{ll} f&{} f\le f_{\text {c}}\, , \\ f_{\text {c}}+k\left[ f-f_{\text {c}}\right] &{} f_{\text {c}}<f\le f_{\text {u}}\, , \\ f^*_{\text {max}}\left[ 1-\exp \left( -a_{\text {f}}\left[ f-b_{\text {f}}\right] \right) \right] &{}f_{\text {u}}< f\, . \end{array}\right. } \end{aligned}$$

(16)

At a critical value $f_{\text {c}}$, the coalescence of microvoids stage is initiated and the void evolution is then controlled by the parameter $k$. Total material failure is described, if the effective void volume fraction $f^*$ attains the value of $f^*_\textrm{f} = 1/q_1$. The proposed modification contains the third case of Eq. (16). Thereby, a robust numerical treatment of total material damage is ensured, for more details see (Seupel et al. 2020; Seupel and Kuna 2019). In this modification, the final effective void volume fraction is set to $f^*_{\text {max}}= 0.995f^*_\textrm{f}$ and the parameters $a_{\text {f}}$ and $b_{\text {f}}$ are prescribed as

$$\begin{aligned} a_{\text {f}}&= \frac{k}{f^*_{\text {max}}-f_{\text {u}}^*},&b_{\text {f}}&=f_{\text {u}}+\frac{1}{a_{\text {f}}} \ln \left( 1-\frac{f_{\text {u}}^*}{f^*_{\text {max}}} \right) \, , \end{aligned}$$

(17)

which ensures a continuous and continuously differentiable transition with respect to $f$. In this formulation, the material is considered as totally failed if the transition value of the effective porosity reaches $f_{\text {u}}^*= 0.98f^*_\textrm{f} < f^*_{\text {max}}$, which corresponds to a transition value of $f_{\text {u}}=\left[ f_{\text {u}}^*+ f_{\text {c}}[\kappa -1]\right] /\kappa$.

2.2 Numerical implementation

The non-local GTN model described in the previous section is implemented in the commercial finite element software Abaqus (2014). Due to the mathematical similarity of the additional Helmholtz-type differential equation (12) for the non-local variable with the stationary heat conduction equation, the implementation can be performed by a User-Defined Material Subroutine (UMAT) (Seupel et al. 2018; Ostwald et al. 2019; Azinpour et al. 2018). This technique avoids the programming of a separate element routine, so that all in-built elements of Abaqus can be used, including features such as contact. The non-local GTN model is therefore implemented in Abaqus/standard as an user defined material UMAT and alternatively for Abaqus/explicit as VUMAT for dynamic problems. In both cases, the constitutive equations are integrated with an Euler backward scheme in the framework of a radial return mapping algorithm including an operator split. The implementation is carried within a monolithic solution scheme for the coupling between nodal displacements and nodal non-local strains. Details about the implementation of the present model can be found in (Seupel et al. 2020). The source code (UMAT/VUMAT) of the present implementation of the non-local GTN model is provided at https://tu-freiberg.de/NonlocalGTN.

3 Simplified parameter identification procedure

3.1 Known trends from literature

A robust strategy for identifying the parameters of the GTN model relies on a profound knowledge of the influence of individual parameters on measurable quantities. For that purpose, the known influences from literature are shortly recalled with respect to the standardized tests to be used in the present approach, i.e., the smooth tensile test and fracture mechanics tests.

3.1.1 Uni-axial tensile test

Regarding tensile tests of common engineering metals like steels, the load-elongation curves are governed by the matrix yield curve ${\bar{\sigma }}\left( {\bar{\varepsilon }}\right)$

$$\begin{aligned} {\bar{\sigma }}&= \sigma _{\text {t}}\approx \sigma \left[ 1+\varepsilon \right] \, ,\ \varepsilon _{\text {t}}=\ln \left( 1+\varepsilon \right) \, , \ {\bar{\varepsilon }}= \varepsilon _{\text {t}}-\frac{\sigma _{\text {t}}}{E} \, . \end{aligned}$$

(18)

The true stress and strain, $\sigma _{\text {t}}$ and $\varepsilon _{\text {t}}$, respectively, are calculated with help of the measured engineering stress $\sigma$, engineering strain $\varepsilon$

$$\begin{aligned} \sigma&=\frac{F}{A_0} \, ,&\varepsilon&=\frac{\Delta l}{l_0} \, , \end{aligned}$$

(19)

and the determined Young’s-modulus E, where F is the force, $A_0$ is the nominal cross section of the tensile specimen, $l_0$ is the length of the measurement domain and $\Delta l$ is the elongation.

The relations for the determination of the yield curve above are valid for nearly incompressible materials, up to the onset of necking of the tensile specimen. Beyond this point, correction methods, e.g., Bridgman correction, are needed to determine the true stress–strain behavior. Alternative concepts of extracting the matrix hardening characteristics can be found in the recent literature (Defaisse et al. 2018; Hao et al. 2019; Tu et al. 2019).

Various empirical hardening functions have been proposed. In guidelines (SINTAP 1999; Zerbst et al. 2007) (SINTAP) and in engineering practice (Sherry et al. 2005), the one-parametric power law

$$\begin{aligned} \sigma _{\text {t}}&= {\left\{ \begin{array}{ll} E \varepsilon _{\text {t}}&{} \text {for } \sigma _{\text {t}}<\sigma _{\text {y}}\, ,\\ \sigma _{\text {y}}\left[ \frac{\varepsilon _{\text {t}}}{\varepsilon _{\text {y}}}\right] ^{1/n} &{} \text {for } \sigma _{\text {t}}\ge \sigma _{\text {y}}\, , \end{array}\right. } \end{aligned}$$

(20)

is established. Therein, $\sigma _{\text {y}}$ and $\varepsilon _{\text {y}}$ denote the characteristic stress and the corresponding strain at initial yielding. The hardening exponent $n \ge 1$ is determined with help of a $\sigma _{\text {t}}/ \sigma _{\text {y}}$ vs. $\varepsilon _{\text {t}}/ \varepsilon _{\text {y}}$ plot in a double-logarithmic diagram as the slope for the domain $\sigma _{\text {t}}/ \sigma _{\text {y}}>1$. The strain hardening law can be obtained by taking into account the additive decomposition of the true strain into an elastic and plastic part i.e. $\varepsilon _{\text {t}}=\varepsilon _{\text {t}}^{\text {el}}+\varepsilon _{\text {t}}^{\text {pl}}$, as well as the constitutive law $\sigma _{\text {t}}=E \varepsilon _{\text {t}}^{\text {el}}$ (implicit formulation)

$$\begin{aligned} {\bar{\sigma }}&=\sigma _{\text {y}}\left[ \frac{ {\bar{\sigma }}}{ \sigma _{\text {y}}}+\frac{E}{ \sigma _{\text {y}}} \, {\bar{\varepsilon }}\right] ^{1/n} \, . \end{aligned}$$

(21)

3.1.2 Fracture mechanics tests

A round-robin (Bernauer and Brocks 2002) revealed, that the simulation results of tensile tests, either smooth or notched, are rather insensitive to the intrinsic length. It was thus proposed to incorporate a fracture test for parameter calibration. Fracture tests of ductile materials are usually evaluated in terms of crack growth resistance curves (R-curves) $J=J\left( \Delta a\right)$. The R-curve is specific to each material if specimens with a sufficiently high level of crack tip constraint are employed, such as compact tension C(T)-specimens or deep-cracked single edged bend SE(B) specimens. Such R-curves can be predicted by the GTN-model for each parameter set. Dimensional considerations show that the predicted R-curves are of the form (Seupel et al. 2020; Hütter et al. 2013a; Xia and Shih 1995)

$$\begin{aligned} \frac{J}{\sigma _{\text {y}}l_{\text {nl}}}&=\mathop {f}\left( \frac{\Delta a}{l_{\text {nl}}},\frac{E}{\sigma _{\text {y}}}, f_0, f_\textrm{n}, f_{\text {c}}, f_{\text {u}}, f_{\text {u}}^*,\varepsilon _\textrm{n},s_\textrm{n}, q_1, q_2,\right. \nonumber \\&\quad \left. \text {dimensionless hardening parameters} \right) \, . \end{aligned}$$

(22)

The fracture toughness $J_\textrm{Ic}$ forms a characteristic point of the R-curve (defined slightly differently in different standards), so that $\Delta a$ drops out from the dependencies of Eq. (22)

$$\begin{aligned} \frac{J_\textrm{Ic}}{\sigma _{\text {y}}l_{\text {nl}}}&=\mathop {f}\left( \frac{E}{\sigma _{\text {y}}}, f_0, f_\textrm{n}, f_{\text {c}}, f_{\text {u}}, f_{\text {u}}^*,\varepsilon _\textrm{n},s_\textrm{n}, q_1, q_2, \right. \nonumber \\&\quad \left. \text {dimensionless hardening parameters} \right) \, . \end{aligned}$$

(23)

And remarkably, the intrinsic length $l_{\text {nl}}$ drops out when considering the dimensionless slope of the R-curve, the so-called tearing modulus,¹ at a characteristic point (Paris et al. 1979)

$$\begin{aligned} T_{\textrm{R}}&:= \left. \frac{E}{\sigma _{\text {y}}^2}\frac{\textrm{d}J}{\textrm{d}{\Delta a}}\right| _{\mathrm c}\nonumber \\&=\mathop {f}\left( \frac{E}{\sigma _{\text {y}}},f_0, f_\textrm{n}, f_{\text {c}}, f_{\text {u}}, f_{\text {u}}^*,\varepsilon _\textrm{n},s_\textrm{n}, q_1, q_2, \right. \nonumber \\&\quad \left. \text {dimensionless hardening parameters} \right) \, . \end{aligned}$$

(24)

Table 1 Fixed ad-hoc GTN parameters

Even in this dimensionless form, these relations comprise nine parameters $\bigl \{f_0$, $f_\textrm{n}$, $f_{\text {c}}$, $f_{\text {u}},$ $f_{\text {u}}^*,$ $\varepsilon _\textrm{n},$ $s_\textrm{n}$, $q_1, q_2 \bigr \}$ which cannot be solely determined from the tensile test. To the authors’ knowledge, there is no study available in the literature, which varied all these parameters simultaneously in a systematic way. Rather, $f_0$ and/or $f_\textrm{n}$ are usually determined from metallographic considerations like Franklin’s formula. The parameters $q_1$, $q_2$ are extracted from micromechanical unit cell calculations (Kuna and Sun 1996; Faleskog et al. 1998; Koplik and Needleman 1988). Mostly the seminal ones of Tvergaard and Needleman are applied ($q_1=1.0$, $q_2=1.5$) as can be seen in the compiled table of employed GTN parameters by Zhu et al. (2022). The nucleation parameters $\varepsilon _\textrm{n}$ and $s_\textrm{n}$ can be determined from metallographic investigations of interrupted tests or partial unloading (Zhang et al. 2021). Often, these values are set ad-hoc to $\varepsilon _\textrm{n}=0.3$ and $s_\textrm{n}=0.1$ (compare also the compilation of parameters in Zhu et al. (2022)). This choice is robust from a numerical point of view and related predictions yielded reasonable agreement with respective experiments. The nucleation porosity $f_\textrm{n}$ is however neglected in this work, as it is not necessarily needed to predict ductile failure at moderate to high triaxialities. Sensitivity studies regarding the R-curve, Eq. (22), were performed by cell model simulations. Gao et al. (1998) and Xia and Shih (1995) varied the initial porosity $f_0$ at a fixed value of $f_{\text {c}}$ and proposed to calibrate $f_0$ to the experimentally measured tearing modulus. Vice versa, a considerable influence of $f_{\text {c}}$ on the tearing modulus at fixed $f_0$ was noted in (Bernauer and Brocks 2002; Hütter et al. 2013a; Chahboub and Szavai 2019; Chen et al. 2020) and it was thus suggested to calibrate $f_{\text {c}}$ to the respective experimental data. Subsequently, the intrinsic length $l_{\text {nl}}$ can be determined from the measured fracture toughness $J_\textrm{Ic}$. Recently, Seupel et al. (2020) found that the same R-curve can be obtained with different sets of GTN-parameters and that it is mainly the ratio $f_{\text {c}}/f_0$ that determines the tearing modulus. A non-uniqueness of the GTN-parameters with respect to macroscopically measurable quantities was also found in (Chahboub and Szavai 2019; Brinnel et al. 2015). Unfortunately, the experimental determination of the tearing modulus is not part of common standardized fracture mechanics procedures and thus not available in most material data sheets. This problem applies to many non-standard tests such as notched tensile tests, interrupted tests or certain quantitative metallographic quantities. Furthermore, all aforementioned calibration procedures require several iterations over numerical simulations to extract relationships like (22)–(24) within an optimization loop.

3.2 Concept of the simplified procedure

The aim of this work is to propose and validate a simple and robust iterative-free identification procedure for the non-local GTN parameters based on a number of certain pragmatic assumptions and the requirements of only two standardized tests. The main assumptions and requirements will be introduced in the next sections, as well as how to determine the required key input parameters and followed by a guideline for the application of the procedure.

3.2.1 Basic assumptions and requirements

Firstly, several parameters are fixed ad-hoc in the usual way as given in Table 1. This pragmatic choice corresponds to typical values, as can be found, e.g., in the collected parameter sets in (Chen et al. 2021; Zhu et al. 2022; Chahboub and Szavai 2019).

As mentioned previously, only two standard tests are required for the completeness of the proposed procedure:

• 1. A uni-axial tensile test (to determine yield stress and hardening exponent n)

• 2. A fracture test C(T) or SE(B) (to determine $f_{\text {c}}$ and $l_{\text {nl}}$)

In this work, the established one-parametric power law (20) is chosen as an appropriate strain hardening law, which is widely accepted within the engineering community (SINTAP 1999; Zerbst et al. 2007) (SINTAP) and (Sherry et al. 2005). The elastic–plastic parameters to be determined are thus reduced to only two, i.e., the power law exponent n and the yield stress $\sigma _{\text {y}}$. For practicality, two possibilities to determine n are presented, see Fig. 1a), depending on the available data from the tensile test:

• (a) Determine n from a logarithmic regression of the experimental stress–strain curve, or

• (b) Estimate n conservatively using $n = \frac{2}{1-R_\mathrm {p0.2}/R_\textrm{m}}$ according to the SINTAP guideline,

where $R_\mathrm {p0.2}$ and $R_\textrm{m}$ correspond to the yield strength (proof strength at $0.2\%$ plastic strain) and ultimate tensile strength, respectively. Consequently, the initial yield stress $\sigma _{\text {y}}$ can then be evaluated from $R_\mathrm {p0.2}$ for the determined n using the following relation

$$\begin{aligned} \sigma _{\text {y}}= \frac{\left[ R_\mathrm {p0.2}\right] ^\frac{n}{n-1}}{\left[ R_\mathrm {p0.2} + 0.002 \, E\right] ^\frac{1}{n-1}} \end{aligned}$$

(25)

Fig. 1

Schematic representation of the necessary data required for the proposed simplified identification procedure

Moreover, according to the dimensional considerations in the previous section, the fracture toughness $J_\textrm{Ic}$ is linearly related to the intrinsic length $l_{\text {nl}}$ but nonlinearly to the critical void volume fraction $f_{\text {c}}$. The latter is related to the tearing modulus $T_{\textrm{R}}$ according to (24), which is independent of $l_{\text {nl}}$. In line with these considerations, two characteristics are required from the R-curves of standard fracture tests (here C(T) following the ASTM E1820 guidelines):

• The tearing modulus $T_{\textrm{R}}$ (to determine $f_{\text {c}})$

• The fracture toughness $J_\textrm{Ic}$ (to determine $l_{\text {nl}}$)

The fracture toughness, according to the ASTM E1820 standard, is determined from the intersection of the 0.2 mm offset blunting line with the power-law fit of the R-curve, see Figs. 1b or 2.² However, the experimental determination of the tearing modulus has not yet been standardized; therefore, it remains an open point to be explored in the next section.

3.2.2 Identification of a robust measure of tearing

Since no suitable measure for the tearing modulus has been given in the literature, two possible choices are proposed and investigated here, in order to find a robust and reliable method of identifying the tearing modulus. Based on the ASTM E1820 standard, the experimental R-curves are usually fitted by a power-law curve of the form $J = C_1 (\Delta a)^{C_2}$, with the dimensionless exponent $C_2$ determining the steepness of the curve, which makes it a possible candidate to quantify tearing. The other feasible choice to define tearing is described by the dimensionless slope of the R-curve as in Eq. (24), taking, however, a linear regression of all experimental points between the 0.2 mm and 1.5 mm exclusion lines, instead of one characteristic point. Figure 2 shows a schematic representation of the two possible candidates for the tearing measures.

Fig. 2

Schematic representation of the two methods to determine the tearing modulus

After introducing the two possible choices for the tearing modulus, it remains to assess which one of them provides a more robust and reliable measure. For that purpose, sets of virtual tests were conducted to show the effects of some of the key parameters on these measures. Based on the dimensional consideration (24), it was found that the tearing modulus should be independent of the intrinsic length $l_{\text {nl}}$. So, as a first test, an arbitrary parameter set is chosen, where all parameters are fixed while the intrinsic length $l_{\text {nl}}$ is varied. Figure 3 displays the R-curves of corresponding virtual C(T) simulations for three different $l_{\text {nl}}= 0.05, 0.1$ and 0.2 mm.

Fig. 3

R-curves for the different intrinsic lengths with the associated different measures of tearing: $C_2$ as the exponent of the power law fit based on the ASTM 1820 guideline and $T_\textrm{R}$ as the linear regression between the exclusion lines

Additionally, the tearing behavior is evaluated using the two different methods described before. From the results in the legend of Fig. 3 it is observed that the tearing modulus $T_\textrm{R}$, if determined by a linear regression line, is indeed independent of the intrinsic length $l_{\text {nl}}$ within acceptable tolerance. In contrast, the R-curve exponent $C_2$ exhibits an inaccaptable scatter of more than 50%. Consequently, the exponent $C_2$ has been dropped as candidate for the tearing measure.

Another set of virtual C(T) tests was conducted to investigate the influence of the hardening exponent n. The results in Fig. 4 show that the tearing modulus $T_\textrm{R}$, again measured via a linear regression, even gave the same results for the different hardening exponents.³ This behavior was confirmed by virtual tests with other values of $E/\sigma _{\text {y}}$ (not shown here). It can be therefore established that the tearing modulus, defined as the linear regression of the R-curve between the exclusion lines, can be considered a suitable and robust measure for tearing. Furthermore, it can be concluded that the tearing modulus can be used to identify the critical porosity $f_{\text {c}}$, thus settling all the necessary requirements for a complete description of the proposed simplified procedure.

Fig. 4

R-curves for the different hardening exponents n with the different measures of tearing

3.2.3 Identification procedure usage guidelines

Based on these findings and the aforementioned literature survey, an iteration-free calibration procedure is proposed which requires only the experimental data of a tensile test and a standardized fracture test. The main steps for the utilization of the procedure are summarized as follows:

• Gather all key parameters: E, $R_\mathrm {p0.2}$, $R_\textrm{m}$, $J_\textrm{Ic}$, and $T_\textrm{R}$

• From $R_\mathrm {p0.2}$ and $R_\textrm{m}$ determine n and $\sigma _{\text {y}}$

• For $T_\textrm{R}$ use the look-up “$T_\textrm{R}$-diagram” to read of $f_{\text {c}}/f_0$ as in Fig. 5a

• For the identified $f_{\text {c}}/f_0$ determine $l_{\text {nl}}$ using $J_\textrm{Ic}$ and the look-up “$l_{\text {nl}}$-diagram” as in Fig. 5b

The whole procedure is also summarized in Fig. 5, which represents a full schematic illustration of the identification procedure. Figure 5a summarizes (in an algorithmic setting) the main ingredients and steps of the identification procedure. A typical example of the look-up diagrams used in this procedure is shown for a specific material strength characterized by the ratio $E/\sigma _{\text {y}}$ in Fig. 5b and c, with the four steps to read off the parameters (in Fig. 5, curves for $E/\sigma _{\text {y}}= 500.0$ are shown as a representative example). A flowchart of the procedure is additionally shown in Fig. 5e, which helps to easily visualize the full parameter identification procedure. The diagrams of Fig. 5 visualize the dimensionless relations (23)–(24) between the GTN parameters and the measurable quantities of fracture initiation toughness $J_\textrm{Ic}$ and tearing modulus $T_\textrm{R}$. Each data point therein was extracted from the simulation of a C(T) test with respective GTN parameters.

Fig. 5

Schematic illustration of the identification procedure, where a summarizes the main ingredients and step of the procedure, b and c represent the typical diagrams used to identify the parameters and d shows the flowchart of the identification procedure

Fig. 6

Example of the FE-mesh of the C(T)-specimen. A radius $r_t$ is applied at the crack tip and a fine elements size $b_e$ is used in the fracture process zone

Figure 6 shows an FE mesh of the standard size C(T)-specimen used for the simulations to create the diagrams. An initial crack length $a_0$ was assumed with $a_0/W = 0.6$, where $W = 50.0$ mm corresponds to the width of the specimen. A 2D plane strain model was employed with an effective thickness of $B_\textrm{eff}= \sqrt{B \, B_\textrm{n}}$. The thickness and net thickness due to side grooves were prescribed as $B=25.0$ mm and $B_\textrm{n}=20.0$ mm, respectively. Due to symmetry, only the half model was used for the simulation. A displacement-controlled loading was applied via an elastic wedge at the pin holes (Linse et al. 2012). A small radius of $r_t = l_{\text {nl}}/2$ was prescribed at the initial crack tip in order to avoid large element distortion (Hütter et al. 2013a). According to previous convergence studies w.r.t the element size $b_e$ (Linse et al. 2012; Seupel et al. 2020), the crack ligament was meshed by square shaped elements (serendipity Q2-Q1 formulation with reduced integration, Abaqus type CPE8RT), with $b_e = l_{\text {nl}}/4$, see magnified view in Fig. 6.

The simulated R-curves were evaluated according to ASTM E1820 (ASTM 2020). Here, the required crack extension $\Delta a$ is defined as

$$\begin{aligned} \Delta a = \Delta a_\textrm{blunt} + \Delta a_\textrm{fail} \,, \end{aligned}$$

(26)

where the first term $\Delta a_\textrm{blunt}$ corresponds to the blunting contribution of the crack tip and is measured as the relative displacement along the ligament between the nodes highlighted by $``\times ''$ in Fig. 6. The second term $\Delta a_\textrm{fail}$ denotes the crack growth due to total material failure along the ligament.

Table 2 Material properties of the three steels under consideration

Fig. 7

Stress–strain curves for the StE 460 steel. Experimental data (Brocks et al. 1995) compared to the power-law fit up to $\varvec{R}_m$ and the power law strain hardening curve w.r.t SINTAP

For each $E/\sigma _{\text {y}}$ ratio, two diagrams are needed, the “$T_\textrm{R}$-diagram” for the tearing modulus $T_\textrm{R}$ (see Fig. 5b) and the “$l_{\text {nl}}$-diagram” for the internal length $l_{\text {nl}}$ (see Fig. 5c). In order to make this procedure valid for a wide range, these diagrams have been created by extensive FE simulations of the C(T) fracture test for many different values of the ratio $E/\sigma _{\text {y}}$ and strain hardening exponent n, to cover the majority of tensile properties of ductile metals used in different applications (The chosen values follow mainly the ranges specified in Sherry et al. (2005)). As has been shown in the previous section, the tearing modulus does not depend on the hardening exponent n, which is also clear from Fig. 5b, where $T_\textrm{R}$ is the same for the different n. So to identify the ratio $f_{\text {c}}/f_0$, one has to draw a horizontal line at the given $T_\textrm{R}$ (step 1 in Fig. 5b), and read directly the ratio $f_{\text {c}}/f_0$ as the intersection point between the horizontal line and the curve (step 2 in Fig. 5b). Next, to identify the internal length $l_{\text {nl}}$, one has to use the the “$l_{\text {nl}}$-diagram” Fig. 5c, by drawing a vertical line at the identified $f_{\text {c}}/f_0$ (step 3 in Fig. 5c), and read the ratio $l_{\text {nl}}\sigma _{\text {y}}/ J_\textrm{Ic}$ as the intersection point between the vertical line and the curve corresponding to the given hardening exponent n (step 4 in Fig. 5c). Knowing now both $\sigma _{\text {y}}$ and $J_\textrm{Ic}$ for the material, the non-local length $l_{\text {nl}}$ follows immediately.

Respective diagrams for other combinations of relative strength $E/\sigma _{\text {y}}$ and hardening exponent n can be found in Appendix A. Therein, it was was taken into account that the tearing modulus does only negligibly depend on n, so that a single diagram, Fig. 15, is required for the relation between $T_\textrm{R}$ and $f_{\text {c}}/f_0$.

4 Benchmark problems

The proposed identification procedure will be tested in this section to evaluate its robustness and reliability. Benchmark problems will be carried out for three different materials, whose respective experimental data are available in published literature. The materials chosen for the benchmark test possess distinct characteristics and are used in different engineering applications: The first material is a mild steel (StE 460) (Brocks et al. 1995), the second one a pressure vessel steel (18Ch2MFA) (Seupel et al. 2020; Abendroth and Kuna 2003; Müller 1999; Abendroth 2005) used mainly in nuclear reactors, the third, a low alloy steel used for forged turbine shafts (27NiCrMoV 15-6) (Kulawinski et al. 2021). The elastic–plastic properties and the fracture parameters for the three materials are summarized in Table 2. Different verification and validation tests will be conducted for each of these materials depending on the available experimental data. Details on the experimental procedures can be found in the given references.

4.1 Mild steel

Brocks et al. (1995) investigated a high-strength fine grained structural steel with old German designation StE 460 (corresponding to S 460N or P 460N in current European standards). Among others, a smooth tensile test, a fracture test with a C(T) specimen and a fracture test under a low level of in-plane constraint with a M(T) specimen were performed. The experimental stress–strain data for this material are available and given in Fig. 7, along with the power law fit and the power law strain hardening curve w.r.t the SINTAP approximation of n. As observed before and as can be seen from the different stress-curves, the SINTAP approximation underestimates the experimental results and gives a rather conservative approach. For all the materials, both methods of estimating the hardening exponent will be employed for the benchmark tests. For simplicity, these two methods will be referred to throughout the remaining paper as follows:

• For the regression of the experimental data points up to $R_\textrm{m}$: “$\varvec{R}_m$-Fit”

• For the SINTAP estimation: “SINTAP”

Table 3 Sets of calibrated parameters for the StE 460

Fig. 8

Predicted C(T) $R-$curves using the parameter sets identified by the proposed procedure compared to the experimental results (Brocks et al. 1995) for the StE 460 steel

For the $\varvec{R}_m$-Fit method, a hardening exponent of $n=4.5$ and a respective initial yield stress of $\sigma _{\text {y}}= 384.0$ MPa were determined, which yield a ratio of $E/\sigma _{\text {y}}\approx 550.0$. Now putting the procedure into action:

• 1. Using the $T_\textrm{R}$-diagram Fig. 15 for the $E/\sigma _{\text {y}}= 550.0$ curve to read the ratio $f_{\text {c}}/f_0$

• 2. Using the $l_{\text {nl}}$-diagrams Fig. 16 for $E/\sigma _{\text {y}}= 550.0$ and $n=4.5$ to determine $l_{\text {nl}}$

The same procedure is performed for the SINTAP method and the identified parameters for both methods are summarized in Table 3.

For the purpose of verification, FE-simulations of the C(T)-specimen using the sets of identified parameters were performed. The predicted R-curves from the simulations were then evaluated following the ASTM E1820 guidelines. The predicted $R-$curves are plotted and compared to the available experimental data in Fig. 8. Both simulated $R-$curves fit the experimental results quite well. It can also be noticed, as expected, that the $R-$curves predicted with the SINTAP method represent a more conservative prediction, characterized by the lower values of the J-integral in the initial stage of the $R-$curve. In this first benchmark test, the identification procedure presents a reliable prediction of the $R-$curves for both methods of determining the hardening exponent, even though the determined n-values for each method lead to different $E/\sigma _{\text {y}}$-ratios, which necessitates performing the identification procedure using different diagrams for each method. The obtained $f_{\text {c}}/f_0$-ratios for the given $T_\textrm{R}$ will therefore not be identical, since the tearing modulus depends on the $E/\sigma _{\text {y}}$-ratio. This actually indicates that despite the differences between the two methods for the determination of the hardening exponent, the proposed identification procedure will still achieve good results, if the right corresponding diagrams are used to read and identify the different parameters.

These results of the C(T)-specimen $R-$curves work very well as a verification test for the proposed procedure, especially since all the diagrams used for identification were generated based on virtual C(T)-tests. Nevertheless, a validation test is still required for a complete proof of concept. For that purpose, the M(T)-specimen was chosen for the StE 460. Following Brocks et al. (1995), a 2D FE model was utilized with plane strain assumption due to the side grooved M(T)-geometry. An FE discretization similar to the requirements of the C(T)-specimen was employed, as previously explained in Sect. 3.2.3. The $R-$curves for the M(T)-specimen were evaluated according to the procedure described in (Aurich 1993), using both identified parameter sets and compared to the available experimental data. Figure 9 shows the different R-curves of the M(T)-test, where the predicted $R-$curves are observed to give a reasonable qualitative and quantitative fit of the experiments.

Fig. 9

Predicted $R-$curves for the M(T)-specimen using the parameter sets identified by the proposed procedure in comparison to the experimental results (Brocks et al. 1995) for the StE 460 steel. C(T) results are plotted to highlight the apparent constraint effects

In addition to the conclusion that the results of the M(T) test validate the concept of the proposed identification procedure, they also highlight another very valuable point. They demonstrate, in fact, that the non-local GTN model is capable of capturing crack-tip constraint effects. The M(T)-specimen is known to exhibit a lower in-plane constraint compared to the C(T)-specimen, which manifested itself in the difference of the steepness of the predicted $R-$curves slopes. The C(T)-specimen, with higher constraints or triaxiality at the tip, corresponds to rather flat $R-$curves compared to those of the M(T)-specimen, which themselves exhibit a much lower triaxiality at the tip. Moreover, all these results underline the practicality of the proposed identification procedure. It should be emphasized that different tests, for one material, can be predicted with reasonable accuracy using the identified parameter set without any additional parameter tuning.

4.2 Reactor pressure vessel steel

For the second benchmark test, the reactor pressure vessel steel 18Ch2MFA investigated by Seupel et al. (2020), Abendroth and Kuna (2003), Abendroth (2005) and Müller (1999) was chosen. The material properties needed for the application of the procedure are summarized in Table 2. The stress–strain curves of the experimental data and the power law curves of the two different hardening exponents are plotted in Fig. 10. For this material, the two hardening exponents lead to a nearly the same $E/\sigma _{\text {y}}$-ratio, so that the same corresponding two diagrams can be used for the identification procedure. Consequently, for both exponents n, the obtained ratio $f_{\text {c}}/f_0$ is the same, since the tearing modulus, as shown before, is independent of the hardening exponent.

Fig. 10

Stress–strain curves for the 18Ch2MFA steel. Experimental data (Seupel et al. 2020) compared to the power-law fit up to $\varvec{R}_m$ and the power law strain hardening curve w.r.t SINTAP

The identification procedure is applied as described before and the resulting parameter sets are given in Table 4. Similar to the previous benchmark test, the C(T)-specimen was used for the verification of the identified parameter sets. The corresponding predicted $R-$curves were evaluated and compared with the experimental results as shown in Fig. 11.

Table 4 Sets of calibrated parameters for the RPV-steel

For this material, both identified sets deliver a very close fit of the experiments. A reason for this accuracy can be attributed to the obtained $E/\sigma _{\text {y}}$-ratio for both methods, which is very close to the given experimental $E/R_\mathrm {p0.2}$ ratio. The $R-$curves predicted with the SINTAP method are again slightly more conservative.

Besides the verification tests, some validation examples were also conducted for this material for the sake of completeness of the benchmark test. For this case, the failure behavior of smooth and notched tensile tests, as well as the Small Punch Test (SPT) were chosen. Details on the underlying FE models can be found in (Seupel et al. 2020). Firstly, the predicted force-diameter reduction response of notched tensile specimens with two different notch radii (R2 and R4) and one smooth tensile specimen are compared to the experimental results in Fig. 12a. Secondly, the force-displacement curves of the SPT are compared to their respective experiments in Fig. 12b.

From the predicted tensile test results in Fig. 12b, it can be seen that the simulation results are in a very good agreement with the experimental results for both the $\varvec{R}_m$-Fit and the SINTAP method. The maximum forces are accurately predicted for all specimen types with the $\varvec{R}_m$-Fit method, whereas, as expected, more conservative values are obtained with the SINTAP method. The failure points are slightly overestimated ; however, they remain in the acceptable range of error.

Fig. 11

Predicted C(T) $R-$curves using the parameter sets identified by the proposed procedure compared to the experimental results for the RPV steel (Müller 1999)

Considering the SPT results in Fig. 12b, the predicted force-displacement curves show a very reasonable match of the experimental data. Moreover, the results with the $\varvec{R}_m$-Fit method are able to capture very closely the maximum force values of the upper bound (UB) of the experiments, and the SINTAP results again yield a conservative predictions. The load drop and failure points match the experimental points quite well for both methods.

Fig. 12

Validation tests for the RPV steel. Experiments by Abendroth and Kuna (2003), Abendroth (2005)

In conclusion, the proposed identification procedure was also able to successfully predict the distinct behaviors of different tests and specimens for this material.

4.3 Forged alloy steel

The final material of choice is the forged-alloy steel 27NiCrMoV 15-6 (also named EN 1.6957) investigated by Kulawinski et al. (2021), which is used for turbine shafts. As for the two previous tests, the hardening exponent n for both methods is first determined and the stress–strain curves of the experiments along with the two power-law curves are plotted in Fig. 13. The material properties are summarized in Table 2 and the identified parameters are given in Table 5.

Fig. 13

Stress–strain curves for the 27NiCrMoV 15-6 steel. Experimental data (Kulawinski et al. 2021) is compared to the power-law fit up to $\varvec{R}_m$ and the power law strain hardening curve w.r.t SINTAP

For this material, the obtained $E/\sigma _{\text {y}}$-ratios were about the same, which means that the same two diagrams can be used for the identification, leading to the same $f_{\text {c}}/f_0$-ratio for both hardening exponents n. The verification tests are once again done with the help of the C(T)-specimen, and the comparison between the experimental results and the predicted $R-$curves are shown in Fig. 14a. A good match between the predicted results and the experiments can be observed.

Table 5 Sets of calibrated parameters for the 27NiCrMoV 15-6

As a validation test for this material, the failure behavior of an SPT is investigated and the simulated force-displacement curves are compared to the experiments of Kulawinski et al. (2021), as seen in Fig. 14b. The geometry and assumptions of the FE model are taken from (Kulawinski et al. 2021), whereas the FE mesh is chosen according to the requirements for the non-local GTN model, as proposed by Seupel et al. (2020). The predicted force-displacement curves of the two methods are in an acceptable agreement with the experimental results, and once again the simulation results enclose the experimental data around the point of maximum force, where the maximum force is rather overestimated for the $\varvec{R}_m$-Fit method and as expected SINTAP method yields a conservative result. A reasonable prediction, however, of the load drops and failure points is attained for both methods.

Fig. 14

Verification in terms of the C(T)-$R-$curves and validation in terms of the SPT force-displacement curves for the low alloy steel 27NiCrMoV 15-6. Experiments of Kulawinski et al. (2021)

5 Recommendations for usage of the simplified parameter identification procedure

In the previous section, the proposed procedure was put to the test for three different materials exhibiting a broad range of properties and features. For all three examples, the procedure exhibited a high level of robustness and reliability. This shows that it is possible to identify the non-local GTN parameters in a practical manner, by using the provided diagrams, without the need for a time consuming and complex iterative calibration scheme. This also indicates that, despite the simplifying assumptions made in the beginning, the proposed procedure is still broadly applicable and sufficiently accurate. Nevertheless, inaccuracies can occur if certain considerations are not taken into account. In what follows, some guidelines are given to ensure the best outcomes of the parameter calibration with this approach:

• For the determination of n, two methods were presented in this work. One systematically overestimated the stress–strain curves and the other was more conservative. So if a conservative approach is required, then the SINTAP method is recommended.

• The ratio $E/\sigma _{\text {y}}$ plays a very important role in this procedure. Because of that, special care should be taken when considering this ratio. After determining the hardening exponent n, the corresponding initial yield stress $\sigma _{\text {y}}$ is evaluated. It is therefore recommended to use the diagrams for this specific $E/\sigma _{\text {y}}$-ratio to read off the different parameters.

• It is possible that the obtained $E/\sigma _{\text {y}}$-ratio falls between two available $E/\sigma _{\text {y}}$-ratios’ diagrams. In this case, it is recommended to perform a linear interpolation between the parameters identified in the two diagrams. However, if conservative results are sought, then reading the parameters from the closest higher $E/\sigma _{\text {y}}$-ratio diagrams is advised.

• Similarly, if the determined n falls between two available n, then a linear interpolation between the two corresponding curves is recommended. And for a conservative result, one should use the curve of the next smaller n.

6 Summary and conclusions

In this paper, a simplified parameter identification procedure for the non-local GTN model was proposed. The GTN model is well-established to predict ductile damage and failure. However, it is not widely used in different engineering applications, mainly to avoid the challenging task of calibrating its parameters. This task is deemed complicated and time costly, since it usually requires numerous FE simulations within an iterative optimization scheme. The aim of this work is therefore to provide a robust, reliable and practical method, which facilitates and hopefully promotes the usage of the GTN model. It is an easy-to-use method, based on simplifying assumptions, such as the on-parametric hardening power-law and some a priori fixed parameters, and does not require any additional FE-simulations.

The proposed strategy is an iteration-free procedure and requires experimental data of only two standardized tests, in which the parameters are read from look-up diagrams, that were created based on systematic studies and produced for a wide range of material properties to cover the majority of ductile metals used in different engineering applications. The main assumptions and requirements of the procedure were introduced and the steps how to determine the key input parameters were explained. Moreover, detailed guidelines for the application of the procedure were provided. For verification/validation purposes, three benchmark tests, for three different materials with distinct features, were carried out to assess the robustness and reliability of the proposed procedure. Finally, some key tips and recommendations were given, of how to obtain the best results from the procedure.

The benchmark tests have shown that the identification procedure exhibits a high level of reliability. For all three materials tested, the simulation results were able to accurately match the experimental data using the identified parameters. Since each of the studied materials possess different characteristics, it is expected that the strategy can be employed for most ductile materials. Moreover, some particular conclusions can be drawn:

• The proposed procedure can be applied for different power-law strain hardening exponents n for the same material and still yields sufficient levels of accuracy. Nevertheless, some exponents n (e.g. SINTAP method) provided more conservative results.

• A robust definition of the tearing measure was found, which is reflected in the linear regression of the R-curves between the 0.2 mm and 1.5 mm exclusion lines.

• The tearing modulus $T_\textrm{R}$ was shown to be independent of the internal length $l_\textrm{nl}$, as well as the hardening exponent n. Additionally, $T_\textrm{R}$ is used to identify the key parameter $f_{\text {c}}/f_0$.

• The diagrams used in the identification procedure showed consistent trends for all material strengths $E/\sigma _\textrm{y}$ and hardening exponents n, making the application of the procedure simple and straight-forward regardless of the particular material properties.

• From the different verification and validation tests, it was concluded that using the identified parameters, the specific behaviors of the different tests can be predicted, without the need of further tuning of the parameters.

• It was pointed out that certain considerations should be taken into account to guarantee the highest possible level of accuracy.

Fig. 15

“$T_\textrm{R}$-diagram” for different $E/\sigma _{\text {y}}$ ratios

As with any newly proposed concept, more benchmark and validation tests are required. The proposed procedure will thus be further tested in the future for different materials and applications.

Look-up diagrams

Fig. 16

“$l_{\text {nl}}$-diagrams” for different $E/\sigma _{\text {y}}$-ratios and hardening exponents n