Application of non-associated flow rule for prediction of nonuniform material flow during deep drawing of tailor welded blanks

Abstract

In order to enhance FE prediction capability, researchers are presently showing interest in applications of non-associated flow rule (NAFR) coupled with Hill48 quadratic function in different sheet metal forming operations. In this work, the concept of NAFR based model was implemented for the first time in FE simulation of deep drawing of DP980-IFHS tailor welded blanks (TWBs) of coated and uncoated sheets. The NAFR model was formulated using two approaches: namely, stress-value based Hill48 as yield function and R-value based Hill48 as the plastic potential function in the first approach, and the vice-versa in the second approach. Also, the classical associated flow rule (AFR) based approach coupled with the anisotropic Hill48 yield model was implemented in the FE simulation for comparison purpose. For improving the prediction accuracy, a mixed hardening equation by combining Voce and Swift hardening law was incorporated as the constitutive equation. It was found that FE simulation implementing the NAFR approach predicted the deep drawing behaviour of the parent materials and TWBs more accurately compared to that of the AFR approach.

Introduction

Over the years, several quadratic/non-quadratic and isotropic/anisotropic yield models were proposed by the various group of researchers for the accurate modelling of material behaviour, including Hill series,^1,2 Hosford series, BBC series, Barlat series and many more ³ with a gradual increase in complexity to include more anisotropic information. Both Yld2000-2D ⁴ and BBC2003 ³ yield models, requiring eight anisotropic data, dominates the recent literature on sheet metal forming due to their accurate prediction of material flow and subsequent failure during various complex sheet forming problems. As an alternative approach, researchers implemented the non-associated flow rule (NAFR), where the yield function and the plastic potential function are not assumed to be identical. Instead of the application of classical associated flow rule (AFR), NAFR is implemented to increase the scope of involving more anisotropic parameters.^5–7 One of the pioneering research on NAFR was published by Stoughton, ⁸ which was based on the Hill48 quadratic yield function, and the proposed approach was able to predict variation in R-values and yield strengths along RD (rolling direction), 45° to RD (diagonal direction, DD) and 90° to RD (transverse direction, TD), simultaneously. It is to be noted that the validity of AFR or NAFR cannot be experimentally established directly, and the fundamental theorems, including uniqueness of stress distribution, stability of plastic flow, and uniqueness of a strain and stress state are mostly used to validate AFR assumption. Stoughton argued the validity of NAFR by proving the same fundamental theorem using the proposed NAFR. Recently, an NAFR model based on Hill48 and Yld2000-2d yield functions was developed to successfully predict the numbers of ears in the case of cylindrical deep drawn cups from aluminium alloy AA2090-T3 (6-ears) and AA5042 (8-ears). ⁹ In this research work, a symmetric stiffness modulus was also developed for the NAFR formulation. Although, most of the studies related to NAFR mainly used R-value based function as yielding criterion and stress-based function as potential function. In the study reported by Kim et al., ¹⁰ stress based Hill48 was used as the yield criterion and the Lankford anisotropy based Hill48 was used as plastic potential function to describe the thinning development in the hole expansion experiment.

However, all the studies related to the application of NAFR based material model were primarily conducted on single monolithic materials. Application of NAFR based material models for tailor welded blanks (TWBs) is not reported in the available open literature. It is useful to mention that the TWBs are one of the important breakthroughs for automotive industries in achieving fuel efficiency by reducing car body weight without compromising passenger safety. ¹¹ TWBs are composed of sheet metals of different grades, sheet gauges, and surface conditions, most preferably by laser welding or friction stir welding, to optimise material usage.^12–14 In general, the component blanks of the TWB and the weld zones are modelled using AFR type material models with mostly Hill48^15–17 and Yld89 models^1,18 as the yield functions. Contrary to the recent trends in sheet metal forming simulation, implementation of more advanced material models (Yld2000-2d, BBC2003, etc.) is relatively rare in the case of TWBs^19,20 as the requirement of the number of material data is increased, and this needs more testing.

With the above-mentioned background, it can be comprehended that the NAFR approach offers an alternative to include a greater number of anisotropic material data even with the quadratic Hill48 yield function. Hence, in the present study, the following major objectives are identified:

• The development of the FE model using the Hill48 based NAFR approach in combination with mixed hardening law (Swift-Voce) for prediction of non-uniform deformation behaviour of DP980-IFHS TWBs under deep drawing conditions.

• Understanding the difference in prediction with the implementation of two types of NAFR approaches: viz. the R-value based Hill48 as yield function along with stress based Hill48 as the plastic potential function in the first case, and stress based Hill48 as the yield function in combination with R-value based Hill48 as the plastic potential function in the second case.

• Validation of the FE predicted results with the experimental results from the deep drawing test of DP980-IFHS TWBs such as cup edge undulation, weldline movement and thickness distribution of the cup wall.

Experimentations

Selection of materials and laser welding

In the present study, two different automotive grade steels with a thickness of 1.2 mm were selected. Among the selected steels, one was dual phase steel, namely DP980 and the other one was interstitial free high strength (IFHS) steel. These two materials have different surface conditions, that is, the DP980 steels were galvanised, and IFHS steels were uncoated sheets. The physical properties, including carbon equivalent, hardness and mechanical properties of the steels, are listed in Table 1. The mechanical properties were evaluated along major material directions, that is, RD, DD and TD and the average values are reported.

Table 1.

General and mechanical properties of the selected steels.

Materials	Carbon equivalent	Hardness (HV)	Yield strength (MPa)	Ultimate tensile strength (MPa)	Uniform elongation (%)
DP980	0.454	340	664	1037.0	8.5
IFHS	0.04	100	163	319.8	25.8

These steels are selected due to their varying range of strength and ductility suitable for various applications, from improving crashworthiness to obtaining ergonomic shape. The DP980-IFHS TWB was fabricated using the laser welding method as it offers advantages like high power density, high aspect ratio, ease of automation etc.^21–23 Keeping the sheets in butt joint configuration and employing a CO₂ laser, the TWBs were welded utilising a power of 3 kW, a scan speed of 3000 mm/min and a shielding gas (Argon) with a flow rate of 20 LPM. The weld orientation was kept perpendicular to the RD.

Weld characterisation

Subsequent to the laser welding, the quality of the weld was assessed in terms of metallographic tests and cross weld tensile testing. Samples taken from TWBs were polished and observed under an optical microscope to confirm full weld penetration. Further, the microhardness tests were conducted across the weld region with 500 g load and 15 s dual time. During the microhardness test, it was observed that the hardness increased in the weld region compared to the corresponding parent metals for both the TWBs. The maximum hardness in the weld region was observed to be 398 HV. SEM micrograph confirmed predominant martensitic phases in that zone, which increased the hardness (Figure 1). The weld region was measured to have a width of 1.5 mm and was termed as a hard zone in the subsequent sections. The hardness values were observed to be reduced in the vicinity of the weld region at the DP980 side, which corroborated well with previous studies. ²⁴ This region is termed as the soft zone (SZ).^25,26 The SZ occurs due to the decomposition of the existing martensite phase in the region due to exposure to the temperature closer to but below the AC1 line. The broken martensitic structure was also observed in the present research (Figure 1). However, the micrograph at the above-right delineates a region from the weld zone and the presence of a certain amount of martensite can be seen there. The particular region is closer to the IFHS side, which is very lean in carbon equivalent, and hence hardness is lower in the region. Finally, uniaxial tensile tests keeping the weld perpendicular to the applied load direction (ASTM E8/E8M-11) were performed to estimate the weld strength. The failure was always observed on the IFHS side (Figure 1(b)). This indicated that the SZ at the DP980 side was not detrimental for the TWB.

Figure 1.

(a) Hardness and SEM microstructure and (b) failed tensile test samples.

Deep drawing experiments

The deep drawing experiments were conducted using a flat bottom cylindrical punch with a diameter of 60 mm and a corner radius of 8 mm. The corresponding die diameter was 64 mm with a profile diameter of 10 mm (Figure 2(a)). The punch and the blank-holder for the deep drawing experiments were attached to the fixtures, which were respectively attached to the primary and secondary rams of a hydraulic press of 100-ton capacity. The die was fixed on a fixture attached to the press bed. The relative positioning of the tooling and the fixture design can be visualised elsewhere. ²⁷ Deep drawing experiments of the parent materials were conducted to estimate the limiting drawing ratio (LDR) of the materials. The maximum blank diameters that were successfully drawn for each material are listed in Table 2, along with corresponding blank holding forces. The LDR evaluations for the TWB keeping the weld line position at the blank centre were also performed, and the results can be seen in Table 2. The drawn cups of parent materials and TWB are depicted in Figure 2(b).

Figure 2.

(a) Deep drawing tooling (dimensions are in mm) with the (b) drawn cups.

Table 2.

Limiting diameter and corresponding blank holding force (BHF).

Blank materials	BHF (kN)	Maximum blankdiameter (mm)
DP980	30.1	132
IFHS	13.53	136
TWB	30.1	107

Material modelling

The formulation of the Hill48 in plane stress condition is shown in equation (1), where G, H, F and N are the anisotropic yield coefficients and ${\sigma }_{\mathrm{xx}}$ , ${\sigma }_{\mathrm{yy}}$ and ${\sigma }_{\mathrm{xy}}$ are components of the Cauchy stress tensor in the orthotropic frame. The four anisotropic coefficients of Hill48 can be calculated either from Lankford anisotropy parameters (equation (2)) or using the normalised yield stresses ratios (equation (3)) at different material orientations, and accordingly, the material model is considered as R-value based Hill48 or stress based Hill48.

$${\overline{\sigma }}^2=H({\sigma }_{\mathrm{xx}}-{\sigma }_{\mathrm{yy}}{)}^2+G{\sigma }_{\mathrm{xx}}^2+F{\sigma }_{\mathrm{yy}}^2+2N{\sigma }_{\mathrm{xy}}^2.$$

(1)

$$\begin{array}{c}\multicolumn{1}{c}{\frac{H}{G}=R_0;\frac{H}{F}=R_{90};\frac{N}{F+G}-\frac{1}{2}=R_{45}}\\ \multicolumn{1}{c}{\mathrm{and}G+H=1}\end{array}$$

(2)

$$\begin{array}{c}\multicolumn{1}{c}{G+H=1;F+H=\frac{{\sigma }_0^2}{{\sigma }_{90}^2};F+G+2N=\frac{4{\sigma }_0^2}{{\sigma }_{45}^2};}\\ \multicolumn{1}{c}{F+G=\frac{{\sigma }_0^2}{{\sigma }_{\mathrm{bb}}^2}}\end{array}$$

(3)

It can be seen from the above-mentioned formulations that the Hill48 material model is enabled only to consider four material anisotropies. Therefore, in this study, a non-associated flow rule (NAFR) was implemented, which enables two different yield functions to be incorporated in FE simulation as yield criterion and plastic potential function (Figure 3). Hence, both the definition of Hill48 yield model (R-value based or stress based) can be used at the same time, as a result of which the number of anisotropic parameters used in the FE modelling can be increased. Thus, if R-value based Hill48 is considered as yield function, then the normalised yield stress ratio along RD and three R-values, namely $R_0$ , $R_{45}$ and $R_{90}$ are considered for determining the Hill48 coefficients. Further, as the stress based Hill48 is used as potential function, the normalised stress ratios along RD, DD, TD and BB are considered to calculate the Hill48 coefficients. Thus, if R-based Hill48 is used as yield function along with stress based Hill48 as potential function, total seven number of anisotropic material parameters are utilised; three R values for R-based Hill48 and four normalised stress values for stress-based Hill48. The same explanation can be extended for another combinational case also. Thus, instead of four, seven material anisotropy related information can be incorporated while performing the FE simulations. In the study, stress based Hill48 as yield function and R-based Hill48 as plastic potential function is denoted as (S + R), whereas R-based Hill48 as yield function and stress based Hill48 as the plastic potential function is denoted as (R + S). In the case of NAFR, the yield criterion ( $\mathrm{\Phi }$ ) can be written as equation (4) in terms of yield function ( $f_y$ ) and equivalent stress considering isotropic hardening ( ${\sigma }_{\mathrm{iso}}(\overline{\epsilon }))$ . Further, the normality flow rule can be expressed in terms of potential function ( $f_p$ ) as equation (5).

$$\mathrm{\Phi }=f_y\left(\mathit{\sigma }\right)-{\sigma }_{\mathrm{iso}}(\overline{\epsilon })$$

(4)

$$d{\mathit{\epsilon }}^P=d\lambda \frac{\partial f_p(\mathit{\sigma })}{\partial \mathit{\sigma }}$$

(5)

Figure 3.

Comparison of AFR with NAFR.

In equation (5), $d\lambda$ is the plastic multiplier factor which in the case of AFR can be proved to be equivalent to effective plastic strain increment $(d{\overline{\epsilon }}^P=\frac{\mathit{\sigma }:d{\mathit{\epsilon }}^P}{f_y\left(\mathit{\sigma }\right)}=\frac{\mathit{\sigma }:d\lambda \frac{\partial f_y\left(\mathit{\sigma }\right)}{\partial \mathit{\sigma }}}{f_y\left(\mathit{\sigma }\right)}=d\lambda )$ using plastic work equivalence principle and Euler’s theorem of homogeneous function. However, it can be shown that, $(d\lambda =\frac{f_y\left(\mathit{\sigma }\right)}{f_p(\mathit{\sigma })}d{\overline{\epsilon }}^P)$ as derived in equation (6) for NAFR.

$$\begin{array}{c}{f}_y\left(\mathit{\sigma }\right)d{\overline{\epsilon }}^P=\mathit{\sigma }:d{\mathit{\epsilon }}^P=\mathit{\sigma }:d\lambda \frac{\partial f_p(\mathit{\sigma })}{\partial \mathit{\sigma }}=d\lambda f_p(\mathit{\sigma })\hfill \\ \mathrm{Hence},d\lambda =\frac{f_y\left(\mathit{\sigma }\right)}{f_p(\mathit{\sigma })}d{\overline{\epsilon }}^P\hfill \end{array}$$

(6)

A schematic diagram is shown in Figure 3 to indicate the difference in AFR and NAFR approach in the stress update algorithm employed in FE simulation.

In the case of AFR, the yield function is used to define the yielding limit of the material, and the same function is implemented to estimate the direction of the plastic strain rate increment, whereas in the case of the NAFR two independent functions are applied to define the material yielding and estimating the plastic strain rate increment.

Uniaxial tensile testings of standard ASTM E8/E8M-11 specimens along the three material directions, namely RD, DD and TD were performed to extract the details related to material anisotropy in terms of yield stress and R-values. Both R-values and yield stress values were assessed at strain values corresponding to 50% of the uniform elongation value of that particular metal. ²⁸ Therefore, the R-values were assessed at approximately 5% strain for DP980 and at 15% strain for IFHS. This is done to average out the evolution in anisotropic properties with an increase in plastic deformation. The evaluated values are listed in Table 3. However, yield stress ( ${\sigma }_{\mathrm{bb}}$ ) in the balanced biaxial stress state was also required for obtaining both yield function and potential function as can be seen in equation (3). Hence, a genetic algorithm-based approach ²⁹ was undertaken to estimate $R_{\mathrm{bb}}$ and ${\sigma }_{\mathrm{bb}}$ values.

Table 3.

Normalised stress ratios and R-values of the selected materials.

Materials	${\mathit{\sigma }}_0/{\mathit{\sigma }}_0$	${\mathit{\sigma }}_{45}/{\mathit{\sigma }}_0$	${\mathit{\sigma }}_{90}/{\mathit{\sigma }}_0$	${\mathit{\sigma }}_{\mathit{bb}}/{\mathit{\sigma }}_0$	${\mathit{R}}_0$	${\mathit{R}}_{45}$	${\mathit{R}}_{90}$	${\mathit{R}}_{\mathit{bb}}$
DP980	1	1.017	1.001	0.937	0.853	0.699	0.725	1.175
IFHS	1	0.99	1.04	1.06	1.037	1.303	1.329	0.780

Hardening equations

Relevant hardening parameters were obtained by fitting corresponding hardening equations to the plastic strain-true stress data obtained from the load-displacement data, which was acquired from uniaxial tensile testing of the material. For TWBs, micro tensile ³⁰ samples along the weld were prepared for evaluating the material properties of the weld. Swift-hardening law (equation (7)), Voce-hardening law (equation (8)) and combined Swift-Voce hardening law ³¹ (equation (9)) were used to fit the plastic strain-true stress data collected from tensile testing results. The parameters for the combined hardening laws for all the parent metals and TWBs are reported in Table 4, along with the corresponding fitness (R²) values. It was observed that the combined hardening laws fitted very well with the experimental data compared with the Swift or Voce hardening laws, as can be seen from the reported R² values in Table 4. The fitted curves for a parent metal (DP980) and a TWB (DP980-IFHS) were depicted in Figure 4 to highlight the nature of the fitted curves. Further, the prediction capability of all the hardening laws, including Swift, Voce and combined Swift-Voce laws, were verified with FE simulation of uniaxial tensile test. The FE tensile testing model was prepared in LSDYNA as per ASTM E8/E8M-11 standard. In the material model, three different flow curve responses were incorporated, and corresponding load-displacement curves were recorded separately for each of the input hardening curves. Finally, for IFHS parent material, the load-displacement curves predicted with the three hardening laws, are compared with the experimental data and shown in Figure 5. It can be observed that the predicted load-displacement curve obtained after incorporating the combined Swift-Voce model agreed well with experimental observation compared to the other two predicted curves. A similar observation was also obtained for DP980 material. Hence, in the present study, combined hardening laws were used in the FE simulations.

$$\overline{\sigma }=K{({\epsilon }_0+{\overline{\epsilon }}_p)}^n$$

(7)

$$\overline{\sigma }=A-B(1-\exp (-C\times {\overline{\epsilon }}_p))$$

(8)

$$\begin{array}{c}\multicolumn{1}{c}{\overline{\sigma }=a\times \left\{K{({\epsilon }_0+{\overline{\epsilon }}_p)}^n\right\}+(1-a)}\\ \multicolumn{1}{c}{\times \{A-B(1-\exp (-C\times {\overline{\epsilon }}_p))\}}\end{array}$$

(9)

Table 4.

Hardening parameters of the materials.

Materials	Evaluated parameters for hardening laws							R2-values
	$\mathit{K}$ (MPa)	${\mathit{\epsilon }}_0$	$\mathit{n}$	$\mathit{A}$ (MPa)	$\mathit{B}$ (MPa)	$\mathit{C}$	a	Combined	Swift	Voce
DP980	1772	0.005	0.17	678.9	−430.9	59.66	0.15	0.994	0.96	0.993
IFHS	591	0.003	0.26	193.4	−248.2	7.78	0.71	0.9998	0.9997	0.999
DP980-IFHS	900	0.005	0.01	547.6	−134.6	51.78	0.27	0.999	0.997	0.998

Figure 4.

The fitted hardening models superimposed with the corresponding experimental data (a) parent materials and (b) TWB.

Figure 5.

Comparison of predicted load-displacement curves from different hardening models with the experimental observation for IFHS sheet metal.

Finite element model

The finite element analysis of the deep drawing processes was performed incorporating a half symmetric tooling model for both parent metals and TWBs employing LSDYNA. The deformable blanks were discretised with a mesh size of 1.2 mm using Belytschko-Tsay shell elements,^19,32,33 keeping five through thickness integration points. ³⁴ The inbuilt re-meshing option in LSDYNA was activated to retain mesh quality during excessive deformation. However, in the soft zone and hard zone, the size of the mesh was approximately $0.5\times 1.2$ due to the smaller width of the region. The detailed mesh sensitive analyses based on which the mesh size of the blanks was selected in this study were based on a previous report by the author. ¹⁹ The friction coefficients between the blank and the tools were reasonably assumed by iteratively comparing simulated load-displacement curves with experimentally obtained results (Figure 6(b)). It was found that the suitable friction coefficients between blank-die and blank-binder were 0.12 and 0.22, respectively, for blank materials DP980 and IFHS. Similarly, the friction coefficients between blank and punch were found to be 0.2 and 0.23, respectively, for DP980 and IFHS. The difference in friction coefficients for each material might be attributed to the galvanised surface of DP980 and the uncoated surface of IFHS. Assigning a uniform friction coefficient for both sides of DP980-IFHS blanks did not result in comparable load-displacement curves. Therefore, separate punch, die and binder segments were designed, so that different friction coefficients can be assigned for each side of the TWB (Figure 6(a)). It was observed that the load-displacement curve corroborated well with experimental results when the friction coefficients for blank-die, blank-binder and blank-punch were kept as 0.14, 0.14 and 0.2 respectively for the DP980 side whereas the friction coefficients for blank-die, blank-binder and blank-punch were kept as 0.23, 0.23 and 0.25 for the IFHS side. The change in friction coefficient at either side might be attributed to the strength difference of the component metals of the TWB for which BHF may not act uniformly on both side.

Figure 6.

(a) Half symmetric FE model and (b) the comparison of load-displacement curves with corresponding punch friction ( ${\mu }_P$ ).

In this study, a Fortran code was developed modifying the example code in LSDYNA theory manual to incorporate Hill48 in the form of NAFR along with the combined hardening model. ³³

Results and discussion

Thickness distributions

The thickness distributions were measured for both the parent and TWB deformed cups to comprehend the deformation behaviour by identifying the regions with localised thinning. The thickness distributions of the deformed cups were measured after taking half section of the full drawn cups. For the parent metal, cup sections were taken along the RD, and for TWBs, sections were taken perpendicular to the weld. Also, the thickness was measured using a pointed anvil micrometre with a least count of 0.001 mm. In the case of parent cups, it was observed that the thickness at the cup bottom is almost the same as the undeformed sheets, and the thickness gradually decreased at the cup corner. ³⁵ The sheet material thickness increases at the flange region with drawing operation due to the development of compressive stress. However, the trend differs in the case of deep drawn TWB cups, where extreme thinning (around 14% of the sheet thickness) is observed at the cup bottom of the IFHS side near the weld. Gradual thickening is observed at the cup wall on both sides of the TWB. Near the weldline, the weaker material is excessively stretched as stronger DP980 pulls the weldline towards itself, and the maximum thinning is observed. The nonuniform nature of deformation for TWB cups is clearly observed from the thickness profile of the cups, which also explains the onset of early thinning and lower LDR for the TWB.

The predicted thickness distribution for DP980 and IFHS parent cups are compared with the experimental results in Figure 7. It can be observed that the NAFR with (S + R) combination and the AFR based FE predicted the thickness distribution more accurately than the NAFR with (R + S) combination for DP980 (Figure 7(a) inset). From the Figure 7(a) inset, it can be observed that NAFR with (R + S) predicted more thickness at the same location of the cup wall compared to the other two approaches. For IFHS, predicted thickness distributions from all the simulation approaches (like NAFR with (R + S) or (S + R) and AFR) does not deviate much from each other (Figure 7(a)). In the case of DP980-IFHS cups, the predicted thickness profile from the NAFR with (S + R) combination is also marginally (with 3%−6% better prediction) closer to the experimental observation for both the sides. Overall, it is observed that the simulations based on NAFR with (S + R) combination provides marginally better prediction of the thickness distribution of deep drawn cups from parent metals and also from TWBs.

Figure 7.

Comparison of thickness distributions (a) parent cups and (b) TWBs.

Cup edge profiles

Plastic anisotropy in the metal results in directional flow of the materials during deep drawing, and it causes the earing in the fully drawn cups. ³⁵ The maximum earing is observed in the case of IFHS cups with a 13% variation in the cup height, while for DP980, a 4% undulation in the cup height is observed (Figure 8(a)). The earing heights can be correlated to the respective planer anisotropy $\left(\Delta R\right)$ parameter of the metals. However, deep drawn cups from TWB show severe undulation in the cup height due to nonuniform deformation in the component (Figure 8(b)). It is observed that IFHS experiences a sudden increase in the cup height than the DP980 side near the weld junction. The yield strength of the DP980 is almost four times higher than that of IFHS. Therefore, DP980 always resists deformation more than IFHS in the case of the TWB. While deep drawing, the stress state in the flange region is tensile along the radial direction and compressive in the circumferential direction. DP980 can withstand more compressive stress and resist the circumferential compressive deformation, and hence more compressive hoop stress is induced in the weaker IFHS material. This excessive compressive stress causes the weaker material to elongate radially outwards near the weldline as the thickening of the material is resisted by the blank holding force. Therefore, the material flow in the outwards direction causes the sudden increase in cup height in the IFHS side of the TWB near the weld. Away from the weld zone, the radial tensile force is more predominant, which draws weaker material into the die cavity causing lower height. On the DP980 side, cup height slightly increases in the region away from the weld due to its resistance to draw into the die cavity. Therefore, severe undulation in the cup height is observed in the case of deep drawing of TWBs of dissimilar material combinations.

Figure 8.

Cup height undulation of (a) parent metals and(b) TWBs (angles with respect to RD).

Figure 8 also compares the experimental cup heights with the corresponding simulation results. It can be observed that the predicted cup height results from the FE simulations are slightly different for DP980. On the one hand, a close inspection revealed that the simulation result from NAFR with (R + S) combination was able to predict the trend in cup height undulation quite properly for DP980 with a uniform deviation from the experimental values. However, the NAFR with (S + R) combination under predicted the cup height, particularly at the TD. On the other hand, the predicted cup height from AFR based simulation is close to experimental observation at RD but largely deviated along TD. In the case of IFHS cups, the predicted cup height undulation from simulations with both the NAFR combinations corroborated well with the experimental observation compared to the results obtained from AFR based simulation. In general, it is reported that FE simulation, incorporating AFR with yield function based on stress values, predicts the height variation between peak and valleys of the earing, ²⁹ and AFR with yield function based on R-values predicts maximum cup heights only. Similarly, in this study, simulations with R-value based AFR predicted only the maximum cup height for DP980 as can be seen in Figure 8. However, incorporation of the NAFR based simulation, particularly NAFR with (R + S) combination, enhances the accuracy of the simulation results in terms of both the height variation and peak heights simultaneously. In the case of the TWBs, the undulation in the cup height is properly predicted, including the sudden change in the cup height in the IFHS side close to the weld by all the FE simulations implemented in this present study. However, close observation reveals that predicted cup height with (R + S) combination is slightly better to match the experimental observation and the improvement in earing prediction was almost 20% greater than the other two approaches.

Weldline movement

Similar to the discussion in the previous section, DP980 resists both radial tension and circumferential compressive stresses, and these two are the two predominant stress states during deep drawing. Hence at the cup bottom, weldline shifts towards DP980 from the initial position. Also, the weaker material is stretched at the cup bottom near the weld. At the cup wall, the weldline shifts towards the weaker material as DP980, being more resistive to the compressive stress, pushes the weldline towards the weaker side. The weldline takes a turn towards the DP980 side at the cup edge. Material of the IFHS side expands towards the outside at the cup edge near the weld, as explained in the previous section. This excess material pushes back the weld towards the DP980 as no material is present at the opposite side to hold out against. So, a curvilinear weld movement is observed in the case of deep drawn cups of TWBs. The weld displacement is more at the cup centre and near the top edge of the cup. In Figure 9, the weldline movement from the simulation with (R + S) based NAFR is correlated well with the experimental results.

Figure 9.

Comparison of the weldline movements for DP980-IFHS TWBs.

Stress-strain analyses of deep drawn TWB-cups

The loci of the Hill48 yield function in plane stress condition were plotted for the two parent metals using both stress value based ( $\mathrm{Hill}{48}_{\mathrm{stress}}$ ) and R-value based ( $\mathrm{Hill}{48}_R$ ) approach. A slight difference in the shape of the yield loci derived, based on stress values and R-values, are observed for DP980, whereas, for IFHS, both the yield functions are almost superimposed on each other. This explains the difference in the prediction of the earing and thickness distribution of the deep drawn parent cups. The earing profile is accurately predicted with R-based Hill48 as yield function and stress based Hill48 as the potential function (R + S combination). Proper representation of the R-values by the yield function and corresponding material flow assessed from stress based plastic potential may be attributed to the proper prediction of earing. In addition, it can also be inferred from Figure 10 that the (R + S) combination predicts more effective plastic strain compared to the other combination. The plastic deformation only occurs when the accumulated equivalent stress, based on the yield function for a given stress state, is more than the effective stress estimated from the hardening rule; otherwise, there will be elastic deformation. Thereafter, the effective plastic strain increment is determined using the plastic potential function. It is to be noted that for the same stress increment potential function will be larger in shape compared to the yield function resulting in more plastic deformation. For the same reason, NAFR with (S + R) combination predicts relatively lower plastic deformation. Hence, the earing profile was properly predicted by the (R + S) combination in NAFR, which also overpredicted the thickness values. For better understanding, the contour of principal stress 1 for DP980 obtained at the same time step was compared in Figure 11. The stress value of a particular node is also highlighted in the figures.

Figure 10.

Hill48 loci for (a) DP980 and (b) IFHS.

Figure 11.

Comparison of principal stress 1 for DP980 material (a) NAFR (S + R) (b) NAFR (R + S) and (c) AFR.

According to the stress update algorithm stress increment is expressed as $d\sigma =C(d\epsilon -d{\epsilon }^p$ ). Here $d\epsilon$ is total strain increment and $d{\epsilon }^p$ is plastic strain increment and $C$ is compliance tensor. Therefore, it can be deduced that the more is the increment in plastic strain, the less would be the stress increment. For NAFR with (R + S) approach, the stress value at the indicated node is 509.495 MPa (Figure 11(b)), and the stress values at the same node for NAFR with (S + R) approach is 517.022 MPa (Figure 11(a)). Hence, it can be inferred that NAFR with (R + S) predicts more plastic strain as the rise in stress is lower at the same node, at the same time step. Hence, NAFR with (R + S) combination predicts more thickening at cup wall (i.e. more deformation). Similarly, it can be inferred that more values of the principal stresses indicate a lower amount of plastic strain increment. Hence, the predicted thickness distribution was lower for NAFR with (S + R) combination, but the values were also close to the experimental observation. The same explanation can be extended for AFR based approach too, where the stress value of the same node point is (∼453 MPa).

For the TWB, the predicted strain states, based on NAFR simulations with (S + R) combination, are compared in Figure 12. The strain values at four nodes were also highlighted for better comprehension. From both the figure, it can be observed that the DP980 side experiences less amount of compressive strain at the cup wall compared to the IFHS side, which results weldline movement and subsequent hump at the IFHS side (arrow marked zone in Figure 12). The hump occurs due to excessive radial deformation of IFHS just near the weld as the material flows in an outward direction to satisfy incompressibility conditions. Excessive near plane strain stretching is also observed at the cup bottom compared to the DP980 side, causing the localised thinning (Figure 7). These observations support the explanations given in the Sections 4.1 and 4.2.

Figure 12.

Comparison of (a) major and (b) minor strain state for NAFR with (S + R). The SZ is muted in the images to indicate the weldline.

A close observation on Figures 7(a) and 8(a) reveals that predicted results from both the NAFR based simulations and the AFR based simulations are relatively closer to one another for IFHS material. The difference in predicted results by the different simulation approaches was clearly visible for DP980, but that was not the case for IFHS. This can be attributed to the similar yield locus of Hill48 obtained from R-value and stress valued based calculation for the IFHS materials (Figure 10(b)). The use of similar yield locus as yield function and plastic potential function makes the FE predicted results similar to the results obtained from AFR based simulation. Hence, NAFR based simulation yields enhancement in prediction accuracy if the R-value and stress value based yield locus defer significantly. It is also to be noted that cup height undulation was closely predicted for IFHS cups compared to DP980 cups. This can be attributed to the intrinsic shape of the yield function. Hill48 being a quadratic yield function, has an intrinsic problem to represent material which shows $\overline{R}$ <1 and ${\sigma }_{\mathrm{BB}}>{\sigma }_0$ . ³ Further, DP steels fall in this category of materials^36–38 but not the IFHS steels. Hence, the prediction accuracy for IFHS materials was better in this study.

Conclusions

The following conclusions can be drawn from the present research work.

1. The non-uniform deformation behaviours in terms of undulated cup height, weldline movement and thickness distribution during deep drawing of TWBs were predicted by the FE model, implementing NAFR based approach reasonably well. Application of NAFR enhanced the prediction capabilities, while avoiding the use of advanced yield formulations.

2. The combined Swift-Voce hardening law had been successfully implemented in the FE model. Also, separate friction coefficients on either side of the TWB in contact with rigid tools were implemented. Hence, it can be concluded that the implementation of Swift-Voce combined hardening law and different coefficients of friction improved the prediction of formability results apart from the use of NAFR in the FE model.

3. The FE prediction accuracy was better with the use of NAFR based material model as compared to that of AFR based approach when the R-based and stress-based formulation lead to different yield loci as observed in the case of DP980 material. If the R-based and stress-based yield functions generate similar yield loci, then the FE predicted results were found to be similar in both NAFR and AFR based approaches, as observed in the case of IFHS material.

4. It was found that the NAFR approach with R-based yield function and stress-based potential function (R + S) yielded better prediction of the undulated cup height (for TWB, improvement in earing prediction was almost 20% greater than other two approaches), whereas the (S + R) combination predicted the thickness profile (improvement in thickness prediction 3%–6% better than other approaches) of the deep drawn cup with better accuracy. Therefore, depending upon the material, the NAFR combination should be chosen.