Prediction of Forming Limit Diagram Using the Marciniak-Kuczynski Method for Ti-6Al-4V Using Different Material Models

The forming limit diagram (FLD) is a widely used tool to assess the formability of a metal sheet [1]. The current study aims to investigate the influence of strain rate, material anisotropy and hardening on the FLD of Ti-6Al-4V predicted by the well-known Marciniak-Kuczynski (M-K) method. The tensile data of quasi-static (8 10-5 s-1), intermediate (0.5 s-1) and dynamic experiments (approximately 1000 s-1) on Ti-6Al-4V sheet are available at three different orientations, with respect to the rolling direction: 0°, 45° and 90°. Different hardening models are taken into account. Also, von Mises and Hill yield criterion are considered. The results show that the influence of the hardening law on FLD is significant. In particular, the most conservative limit strains are predicted by the Voce law because of its saturation characteristic. The yield criterion is found to only affect the right part of the FLD. Regarding the strain rate influence, the left part of the FLD is mainly dominated by the amount of uniform elongation, while the right part is strongly dependent on the yield function used. Therefore, for this region the effects of strain rate and yield function are difficult to distinguish. Finally, the effect of material anisotropy on the FLD is significant. Under quasi-static conditions, the Lankford coefficient seems to be the driving factor in uniaxial and equibiaxial deformation. However, in plane strain conditions the effect of the strain hardening exponent is dominant.


Introduction
The assessment of material formability is a crucial requirement in any sheet metal operations. For this purpose, the forming limit diagram (FLD) is nowadays a widely used tool. The FLD is a limit curve separating the principal in-plane strain space into two regions. The region above the curve refers to strain states in which localized necking is expected to occur, while the region below represents the safe region in which the sheet can be formed without failure. Even though the FLD can be experimentally determined, the procedure involves many different specimen geometries and is time consuming. To overcome this issue, several theoretical techniques have been developed over the years. One of the most commonly used was developed by Marciniak and Kuczynski in 1967, namely the M-K method [1].
Since the formability of metals is affected by the imposed strain rate, also the FLD is dependent on it. The FLDs obtained through the M-K technique are able to capture that strain rate sensitivity [2,3]. Also, the influence of anisotropy, introduced in the metal sheet during the rolling process, has been found to be relevant [4]. The M-K method can account for anisotropy, too. Because of its great flexibility, the M-K method has been used for a broad class of materials. In particular, several studies have been carried out on Ti-6Al-4V, the most used titanium alloy. It has been observed that at high temperatures, the effect of the yield function on the FLD is more pronounced than the effect of the hardening model [5]. Other studies highlighted the negative effect of strain rate on the FLD limit strains [6]. However, the anisotropy of Ti-6Al-4V was not taken into account in this work. In general, another important feature of almost every study is the lack of orientation-specific FLDs, namely FLDs along different orientations of the principal strains. The main reason is that the effect of in-plane anisotropy on the FLD is usually neglected in engineering applications, and the FLD is determined with the major principal strain along the rolling direction. The other reason is that the limit strains obtained along different orientations usually lie within the experimental scatter of the FLD, thus the anisotropy influence cannot be clearly distinguished. However, thanks to innovative strain measurement techniques, such as Digital Image Correlation (DIC), it is possible to reduce this scatter and make the anisotropy effect on the FLD relevant [7].
The current work aims to investigate the influence of the aforementioned factors on the predicted FLD of Ti-6Al-4V. The tensile data of quasi-static (8·10 -5 s -1 ), intermediate (0.5 s -1 ) and dynamic experiments (approximately 1000 s -1 ) on Ti-6Al-4V sheet are available at three different orientations, i.e., 0°, 45° and 90°, with respect to the rolling direction. For each of the three strain rate conditions, the FLD is evaluated through the M-K method implemented in MATLAB, considering different material models. More specifically, an isotropic and anisotropic response of the material is assumed using the von Mises and Hill yield functions, respectively. For the latter, two parameters calibration techniques, based on yield stresses and Lankford coefficients, respectively, are considered. The influence of the adopted hardening model is also studied and orientation-specific FLDs are presented. The different FLDs obtained are then compared and critically discussed.

Material and Methods
Material. The material considered in this study is titanium alloy Ti-6Al-4V. It comes in the shape of a cold rolled sheet with a thickness of 0.61 mm. The nominal composition in weight percent is reported in Table 1.  Fig. 1 are obtained for each strain rate and each orientation.

886
Achievements and Trends in Material Forming For each test, the DIC data additionally allowed to derive the Lankford coefficient, according to Eq. 1. εẏy p is the plastic strain rate along the transverse to the loading direction and εżz p the plastic strain rate along the thickness. Thus, the Lankford coefficient is an indication of both the anisotropy and formability of the material. The plastic strain rates of Eq. 1 are obtained considering the slope of the line fitting the plastic strain vs time experimental data in the uniform plastic deformation region. For this reason, the levels of strains between which the Lankford coefficients are obtained depend on the specific plastic curve considered. (1) Material model. For each of the three different strain rates, the plastic curve is fitted through a hardening model for each orientation. Swift, Voce and Combined Swift-Voce hardening laws are considered and reported in Eq. 2, Eq. 3 and Eq. 4, respectively.
contribution of Swift and Voce laws. The Combined model turns out to be the most accurate one in describing any of the plastic curves, thanks to the fitting flexibility provided by its six parameters. Since for each different angle to the rolling direction the fitting parameters are obtained considering the experimental data of that specific orientation, a big amount of parameters sets is available. Thus, in Table 2 only the fitting parameters obtained for the quasi-static experiment along the rolling direction, corresponding to the reference, are presented.
Because the hardening models involved are not strain rate and temperature-dependent, the adiabatic heating, especially occurring at high strain rates, is not explicitly accounted for in the M-K analysis. However, strain rate and thermal softening effects are indirectly considered by fitting the hardening parameters using high strain rate test data. Then, an associated flow rule (AFR), in which the plastic potential coincides with the yield function, is assumed. The governing equation of AFR is reported in Eq. 5, where εṗ is the plastic strain rate tensor, εṗ the plastic equivalent strain rate, f the yield function and σ the stress tensor: In order to take into account the effect of anisotropy, the Hill48 equivalent stress σeq is considered as yield function f and reported in Eq. 6, expressed in the reference frame associated with the axes of orthotropy: F, G, H, L, M and N are the Hill parameters, which are determined, for each strain rate, following two different calibration approaches. The first, here named Sr-based approach, makes use of both yield stresses and Lankford coefficients. The second, here named r-based approach, involves the Lankford coefficients only. The Sr-based and r-based parameter calibration approaches are presented in Table 3, where the subscript of yield stresses σ and Lankford coefficients r refers to the angle with respect to the rolling direction along which they are evaluated. Since this study assumes that the sheet is subjected to plane stress conditions, σyz and σxz are equal to zero. For this reason, the L and M parameters are not involved and their expression thus not reported.

888
Achievements and Trends in Material Forming Hill parameter Sr-based r-based In the description of the isotropic behaviour, the yield function f is the von Mises equivalent stress. It can be obtained by setting F=G=H=0.5 and L=M=N=1.5 in Eq. 6.

Marciniak-Kuczynski method.
Given the hardening law and the flow rule, the formability of Ti-6Al-4V is evaluated by means of the Marciniak-Kuczynski (M-K) method. It assumes an inhomogeneity of the material in the form of an imperfection zone. As shown in Fig. 2, this region, represented by the groove b, is characterized by a thickness lower than the one of the homogeneous zone a. The nature of the imperfection can lie in different factors, such as surface roughness or a macroscopic defect. When subjected to in-plane deformation, zone a deforms proportionally, while in the groove plastic strain accumulates faster because of its thickness reduction. This leads to higher strain increments in the groove, which eventually results in strain localization. The onset of strain localisation is assumed to correspond to the forming limit of the sheet. In the current study, the M-K method is implemented in a MATLAB script, following the iterative algorithm of Fig. 3. Two basic boundary conditions are the fundaments of the algorithm, namely the force equilibrium along the n-axis and congruence between zone a and zone b along the t-axis. Let XYZ and ntZ be the reference frames associated with the loading direction and the groove orientation, respectively. Also, let ψ be the orientation of the groove, namely the angle between the aforementioned reference frames. In the algorithm, firstly the strain path ρ and the groove orientation ψ are set. The former is defined in Eq. 7, where εxx a and εyy a are the principal strains of zone a along the X and Y-direction, respectively.
Then, the strain increment Δεxx a is set, εxx a and εyy a are updated and the strain εzz a is evaluated assuming volume conservation. Therefore, ε a XYZ, the strain tensor of zone a in the XYZ reference Key Engineering Materials Vol. 926 frame, is built. By means of ε a XYZ, the strain tensors ε a ntZ and ε a RTZ are obtained with a change of reference frame. The former is expressed in ntZ reference of frame, the latter in the reference frame associated with the rolling and transverse direction (RTZ). In the current study, it is assumed that RTZ is inclined at an angle α with respect to XYZ. The stress tensor σ a RTZ is then calculated by applying one of the hardening laws of Eqs. 2÷4 and the flow rule of Eq. 5. Afterwards, the congruence condition is applied and the strain εtt b of the groove along the t-direction is obtained. The strains εnn b and εnt b are first guessed and the strain tensor ε b ntZ is then evaluated. By changing the reference frame and applying both the hardening law and the flow rule, the strain tensor ε b RTZ and the stress tensors σ b RTZ and σ b ntZ are obtained. The updated values of εnn b and εnt b are then calculated through the Newton-Raphson method applied to the equilibrium equations. The whole process just described continues as long as the ratio between principal strains in zone b and zone a is lower than a critical value Rlim, here assumed equal to 4. Once Rlim is reached, the strains εxx a and εyy a are considered as the limit values for the specific orientation of the groove.
The iterative procedure of Fig. 3 is repeated, with the same strain path ρ, for different groove orientations ψ and the point reported in the FLD is chosen as the one exhibiting the minimum values of εxx a and εyy a . Finally, if this procedure is repeated as many times as the different ρ values, the FLD is built. Regarding f0, the ratio between the initial thickness of zone b and zone a, it is usually calibrated in order to match the estimated and the experimental values under plane strain conditions [8,9]. However, since these experimental data are not available for Ti-6Al-4V, it is set equal to 0.99.

Results and Discussion
The Lankford coefficients obtained from the quasi-static experiments are represented in Fig. 4a, together with the fitted curve, as a function of the loading angle θ to the rolling direction. The curve is calculated according to the von Mises and Hill yield criterion, the parameters of which are identified using the two different approaches of Table 3. The von Mises curve predicts Lankford coefficients equal to 1 for all the loading directions because anisotropy is not involved in the corresponding criterion. The Hill curves clearly predict different values of the Lankford coefficients. Only the curve evaluated with the r-based approach is able to capture the experimental data. On the other hand, only the Sr-based curve, namely the curve obtained with the Hill function calibrated with the Sr-approach, can predict the yield stresses. The reason why the Hill yield criterion is not able to simultaneously reproduce the experimental behaviour of Lankford coefficients and yield stresses, lies in the nature of the function itself. Because it is quadratic in stresses and orthotropic, a specific relationship exists between yield stresses and Lankford coefficients. From Table 3, for example, if the two different approaches for the definition of the F parameter are combined, the following is obtained: From Eq. 8, it is clear how the Hill criterion is not able to fully predict, for instance, the quasistatic experimental response of Ti-6Al-4V, which exhibits σ0 ~ σ90 but r0 ≠ r90.

. Experimental Lankford coefficients as functions of the loading direction, together with the von Mises and Hill fitted curves, for all the experiments (b) and quasi-static only (a).
In Fig. 4b the same graph of Fig. 4a is shown, but for different strain rates and with the r-based Hill fitted curves only. The quasi-static and intermediate curves exhibit similar values for loading angles close to the rolling and transverse direction. However, they tend to diverge while moving towards θ=45°, where the intermediate curve is characterized by higher values. On the other hand, the dynamic curve shows a significantly different behaviour in the right part of the graph, where it exhibits Lankford coefficient values much lower than the quasi-static and intermediate ones.
By means of the M-K method previously presented, the FLDs can be evaluated considering different strain rates and loading orientations with respect to the rolling direction. Additionally, the influence of the plasticity model, i.e., yield function and hardening model, can be assessed.
For the FLDs of Fig. 5a, the quasi-static hardening properties along the rolling direction are used for the fitting of the hardening laws. Also, the Sr-Hill yield function is considered. The limit strain for plane strain conditions (ρ=0) is quasi-independent from the hardening law. However, while moving towards uniaxial (ρ=-0.5) and equibiaxial deformation (ρ=1), the FLD evaluated with the Voce law tends to diverge from the ones obtained with Swift and Combined laws. The former predicts lower limit strains. The origin of the trend lies in the mathematical expression of the Voce law, see Eq. 3. In fact, the flow stress given by the Voce law tends to saturate for sufficiently high plastic strains, thus predicting larger strains at the same the same stress level. Hence, the FLD is negatively affected. On the other hand, the FLDs obtained with the Swift and Combined law exhibit very similar trends and values, because the corresponding hardening curves nearly overlap.
If the three considered hardening models give a significantly different outcome, as is the case for the 45° direction, the corresponding FLDs themselves can be clearly distinguished, as shown in Fig.  5b. Fig. 5b reports the FLDs obtained at the same strain rate condition of Fig. 5a but using the experiment along the 45° direction to calibrate the hardening laws. In this case, the discrepancy between the limit strains predicted by the Swift and Combined law is much higher than between the corresponding curves of Fig. 5a, due to the pronounced differences in the corresponding hardening curves.
Finally, even though the Voce law predicts more conservative FLDs, the most reliable ones are obtained with the Combined law. Indeed, it is the most accurate in describing the plastic behaviour of the material, thanks to the higher flexibility provided by the six fitting parameters, against the three parameters for both the Swift and Voce laws. For this reason, all the FLDs presented from now on are evaluated using the Combined law. The FLDs obtained using the Combined hardening rule and different yield functions are shown in Fig. 6a, considering the quasi-static experiment in the rolling direction to fit the hardening law.
The left part of the FLD (εmin<0) is not affected by the yield function, as well as the point corresponding to plane strain condition. On the other hand, the right part of the FLDs (εmin>0) exhibits significant differences. Indeed, starting from εmin ~ 0.1, the von Mises criterion leads to limit strains higher than the ones predicted by the Hill yield functions, calibrated with either the Sr-based or rbased approach. The greatest discrepancy between the three models occurs for equibiaxial deformation, where the limit strain predicted by the von Mises yield function is approximately twice the one predicted by the r-Hill function. Moreover, the FLD evaluated with the Sr-Hill function lies in between the two other FLDs in the right part of the graph. This may occur because in quasi-static conditions the anisotropy in stresses, on which Sr-Hill function is based, is less pronounced than the anisotropy in Lankford coefficients, characteristic of the r-Hill function [10].

892
Achievements and Trends in Material Forming Similar trends are also observed in the FLDs based on the intermediate strain rate experiments, here not reported. However, the dynamic FLDs are clearly different, as shown in Fig. 6b. In fact, the Sr-Hill yield function predicts limit strains very close to the ones predicted by the von Mises yield function. This has to be attributed to the anisotropy in stresses observed for the dynamic experiments that is less pronounced compared to the one observed for the quasi-static and intermediate tests.
As previously discussed, the strain rate can have a significant influence on the formability of the material. This is clear in Fig. 7, where the FLDs evaluated for different strain rates and yield functions are shown. Irrespective of the yield criterion involved, the left part of the FLD shows a clear trend. Regarding the right part of the FLD, the considered yield function plays an important role. The formability prediction using the von Mises yield criterion (Fig. 7a) is consistent with previous studies carried out with the same criterion [6]. For the different strain rates, the FLD slope is similar, except for plane strain conditions. Moreover, the quasi-static limit strains are slightly higher than the dynamic ones, which are in turn greater than the ones obtained for the intermediate case. Different strain hardening exponents lie at origin of the differences between the curves. Indeed, several studies have highlighted that high strain hardening exponents delays the onset of localized necking, thus increasing the limit strains [11]. In fact, the highest values of strain hardening exponent and limit strains are observed for the quasi-static experiments, the lowest for the intermediate tests.
If, on the other hand, the Sr-Hill yield function is used (Fig. 7b), the main difference lies in the dynamic FLD, exhibiting the highest limit strains for most of the right part of the FLD. If the r-Hill yield function is used, the dynamic limit strains lie in between the quasi-static and intermediate FLDs (Fig. 7c).  To assess the effect of the anisotropy on the material formability, three orientation-specific FLDs, evaluated under quasi-static conditions, are represented in black in Fig. 8. For each of them, the Combined hardening law corresponding to that specific orientation is considered and the r-Hill yield function is used. Then, the different FLDs are fitted through a MATLAB in-built three-dimensional curve, making sure that all the evaluated points lie on the curve. In the neighbourhood of uniaxial and equibiaxial deformation, the limit strains evaluated for the experiment along the 45° direction are higher than the ones evaluated along the rolling (0°) and transverse (90°) direction. This can be linked with the higher value of the Lankford coefficient observed for the 45° direction experiment. However, close to plane strain straining conditions, the trend is different. Indeed, here the highest limit strain is observed for the FLD along the transverse direction, for which the strain hardening exponent is higher. Thus, the Lankford coefficient seems to control the FLD except for the strain states close to plane strain, where the strain hardening exponent plays the most relevant role.

Conclusions
The forming limit diagrams of Ti-6Al-4V sheet are evaluated using the well-known M-K method. The experimental data used in the analysis come from tensile tests in quasi-static, intermediate and dynamic conditions. Three hardening models are used, namely Swift, Voce and Combined Swift-Voce. An AFR is implemented, making use of von Mises and Hill yield criteria.
The influence of the hardening law on the obtained FLDs is not negligible. In particular, the limit strains predicted by Voce law are always the lowest due to the saturation of the stress inherent to the equation. Despite of the conservative nature of the FLD predicted by Voce hardening, the most reliable limit strains are evaluated by the Combined law, because of its ability to capture the experimental data with high accuracy.
The yield criterion is found to strongly affect only the right part of the FLD. In quasi-static and intermediate conditions, the strains predicted with the von Mises yield function are the highest, with the greatest discrepancy occurring for equibiaxial deformation. On the other hand, in dynamic conditions, the FLDs obtained through von Mises and Sr-Hill yield functions almost overlap, except for the equibiaxial condition. This is because in dynamic experiments, the anisotropy in stresses is less pronounced than in quasi-static and intermediate conditions.
Regarding the influence of strain rate, the left part of the FLD is mainly dominated by the amount of uniform elongation. The higher it is, the higher the limit strains are. However, the trend observed in the right part of the FLD is strongly dependent on the yield function used. Therefore, for this part the effects of strain rate and yield function are difficult to distinguish.
Finally, also the anisotropy of the material affects the FLD. In quasi-static conditions, the limit strains predicted for the 45° direction are the highest in uniaxial and equibiaxial deformation. This may be linked with the highest values of Lankford coefficients observed along that direction. However, the limit strains for plane strain conditions follow the trend of the strain hardening exponent.
To better take into account the effect of anisotropy, more accurate yield functions, able to simultaneously capture the stress and plastic anisotropy, are needed [12]. Moreover, a validation campaign is necessary to further assess the accuracy of the present study.
Key Engineering Materials Vol. 926