Paper received: 17.07.2008 Paper accepted: 25.09.2008 Usage of the Yield Curve in Numerical Simulations Pino Koc1* - Boris Stok2 1 Faculty of Mathematics and Physics, University of Ljubljana 2 Faculty of Mechanical Engineering, University of Ljubljana A comparison of two approaches used in the yield curve characterization of the same material is given in the paper. The first approach is commonly used Ludwig's law with the extension over large strains based on the pre-necking response of a tensile test specimen, whereas the second approach is inverse identification which is based on the post-necking behaviour of the same tensile test specimen. Features of both approaches are examined in the tensile test and deep drawing simulations. In the tensile test simulation the inverse identification method proved to be superior over Ludwig's law. The deep drawing simulation demonstrates how inappropriate yield curve usage leads to wrong predictions. © 2008 Journal of Mechanical Engineering. All rights reserved. Keywords: tensile testing, yield curve, inverse identification metod, deep drawing 0 INTRODUCTION In the contemporary mass production of industrial goods numerical simulations have become an indispensable tool by which product development cost might be substantially lowered. Numerical simulations in a design phase of product development might reduce a number, or even suppress, unsuccessful prototype trials, which in turn can accelerate product launching time. To support the above statements a reliable numerical model is needed which, in general, is not easy to acquire. The so called "default" values of numerical model parameters offered by commercial computer codes are usually not appropriate, especially when modelling high-demanding metal forming processes. Such an easy-to-handle (a credulous) approach might lead to unrealistic results of a performed numerical simulation. To prevent such an undesirable event a great care is necessary when choosing each parameter of the numerical model. An integrated knowledge from various disciplines, such as mechanics, numerics, engineering, technology, material science, etc. is desired to cope successfully with this task. In the paper a small but important fragment of numerical modelling, related in particular to computer simulations of metal forming processes, is addressed: the yield curve usage. In technology processes, where plasticity of the material is of the prime concern, the yield curve is namely a basic material datum that enters in numerical simulations. Consequently, a reliability of those simulations depends significantly on the physical objectivity of the adopted yield curve. The yield curve is obtained by different experiments among which the tensile test is the most widespread in test laboratories. The standardized procedure of the tensile test [1] assumes homogeneous strain and stress fields in the parallel central part of the cylindrical or flat tensile specimen. The logarithmic strain in the longitudinal direction p and true (Cauchy) stress a are obtained from the following equations: i L p = ln— . F A (1a) (1b) where L0 and L represent respectively the initial and the instantaneous length of the inspected domain, measured by an extensometer, F the applied force and A the actual cross-sectional area of the specimen. Then the yield curve can be drawn in a diagram as a plot of a against p. This experimental-analytical method of the yield curve identification is simple and works well as long as the assumption of the stress and strain field homogeneity in the observed domain of the specimen is valid. The yield curve defined up to the onset of necking in accordance with Eqs. (1) is quite sufficient for a great majority of components for which only a mild plastic strain, considering that the corresponding equivalent uniaxial stress-strain states in those structural parts are far below the Corr. Author's Address: University of Ljubljana, Faculty of Mathematics and Physics, Lepi pot 11, SI-1000 Ljubljana, Slovenia, pino.koc@fmf.uni-lj.si 0 a = necking stress and strain, is developed during component forming. But, for some highly strained components, such as deep drawn parts, automotive stamped parts, etc., the yield curve defined in the above manner is not representative after the onset of necking, i.e. after the maximum tensile force F=Fmax is registered. Because the onset of necking gradually changes uniaxial stress state to complex triaxial one, it is physically inappropriate to determine the yield curve by the standard analytical procedure (Eqs. 1). In order to characterize objectively the complete yield curve and not only the portion of the curve prior to the onset of necking an improved description of the material behaviour is needed. In the past many attempts were made to cope with the addressed problem. Mainly, it is derived from the yield curve that is obtained upon the experimental data registered in the tensile test prior to the occurrence of necking by its analytical continuation beyond the necking strains [2] to [4]. Among many commonly accepted functional laws, Ludwig's law is probably the most frequently used nonlinear material law in the engineering practice: a = a + Cq>" (2) In Eq. 2 ay, C and n represent material parameters, which are to be adequately tuned in order to obtain the best fit between a yield curve, which is analytically approximated in accordance with the assumed law, and the pre-necking measured data. Usually, the best fit is obtained by the well known least squares method (LSM). Although no physical coverage in the range of large strains can be attributed to such approach, it is frequently used in computer codes [5] and [6]. It is worth mentioning that the ultimate strain, before tearing of the material commences, can not be known by using this method. Since the termination point is not known, one can extend the yield curve beyond all reasonable limits which can lead, when used in numerical simulations, to unreliable, even unsafe results. To circumvent the described shortcoming of the considered functional extension, additional experimental information regarding post-necking deformation is inevitably to be taken into account. Another approach, also frequently used for estimation of the true stress-strain relation after tensile necking occurrence, is Bridgman's method [7]. It represents an analytical approach based on three assumptions regarding uniform strain distribution, constant principal stress ratio and measurable curvature of the necking region. For the tensile specimen with a circular cross-section the true stress beyond necking is: 1 1 + » 1 ln f 1 + a a F ^min (3) where a and Amin are respectively the radius of the smallest cross-section in the neck and its appertaining cross-sectional area; and R is the radius of curvature of the neck. Bridgman's correction method can calculate the true stress fairly well for a circular cross-section of the tensile specimen. A drawback of the method is its impractical application, since it requires frequent interruption of the tensile test owing to a and R measurements. Based on the same three assumptions Bridgman extended his method to specimens with a rectangular cross-section. However, necking of flat specimens is much more complicated than that of specimens with circular cross-section. Thus the abovementioned assumptions were found to be incorrect [8] and the predicted stress far from the true stress. For this reason the Bridgman method is not used in sheet metal yield curve identification. Alternative to the abovementioned methods are mixed experimental-numerical methods [9], where by means of numerical simulation of a real experiment, the yield curve or a part of it can be identified. Since simulations based on numerical modelling can cope with complicated material laws, complex boundary conditions, etc., experiments are no more limited to the simplest one, e.g. the tensile test prior to necking. Therefore non-trivial, i.e. inhomogeneous, transient and multiaxial mechanical states are allowed in a specimen. By reliable numerical simulation of such an experiment and with appropriate tuning of material parameters in the numerical model, matching between the response of the real experiment and results of the corresponding numerical simulation can be obtained. One drawback of the method is additional software analysis which usually demands adequate mastery. Usage of this method, which is called also the inverse identification method, is nowadays in a constant growing rate and many a = material properties identification problems are solved with it [4] and [9]. In the sequel, two yield curves, obtained from the same tensile test measurements, but each employing a different method of the yield curve definition (functional law extension and inverse identification) will be presented. Then, the curves will be used in two comparative numerical simulations of the same deep drawing process and finally, differences between results of the two simulations will be discussed. 1 YIELD CURVE CHARACTERIZATION In the experimental work, performed in the laboratory for measurement of the Kovinoplastika Lož company, the standard tensile tests of sheet metal were done [1] (Fig. 5). Due to the expected anisotropy, three sets of specimens were machined, each set having a different orientation with respect to the rolling direction of the sheet. The specimens of thickness of 1.22 mm were made of low cost low alloyed steel sheet, used in automotive industry. Results from the tests were processed in a classical manner, i.e. different statistical parameters regarding repeatability and uncertainty were checked, the yield curves according to Eqs. (1) were plotted and orthotropic material parameters were extracted. Fig. 1 represents the yield curves for specimens, whose axis is aligned with the rolling direction of the sheet. One can see that a small difference between the yield curves is found in the mild strain range (say below 9 « 0.25) and that a significant difference appears in the range of strains beyond 9 « 0.25. This difference is 550 11 10 ll F max analytical method inverse identification 0.00 0.05 0.10 0.15 0.20 0.25 0.30 P (/) Fig. 1. Analytically derived yield curves -lengthwise direction 0 5 10 15 20 25 30 ALm (mm) Fig. 2. Measured force of the representative specimen attributed to the material intrinsic damage process, which is highly nonlinear and sensitive to small differences of the material properties, surface quality, etc. of each specimen. The inverse identification procedure, which will be used in the sequel, requires as an input the real measurement data, thus no preprocessing of data in a form of averaging of several measurements is needed. Therefore, data from one single tensile test, which seems to be a rough average of all tests from a particular specimens set, are chosen to be the representative data used for the yield curve identification. Those representative data are shown in Fig. 2 as a plot of the measured tensile force Fm against the extensometer's extension Lm, and in Fig. 3 as a corresponding plot of calculated stress a against strain p (Eqs. 1). The maximum measured force (Fig. 2) has a special importance, because the 550 - 500 ^450 CL 5 k 400 350 300 0 0.05 0.1 0.15 0.2 0.25 0.3 p (/) Fig. 3. Analytically derived yield curve of the representative specimen 9 8 7 onset of necking is usually associated with it. Also, with regard to thus obtained stress-strain relationship this value sets the limit on the yield curve, beyond which validity of stress is lost (Fig. 3). The applicable part of the yield curve in Fig. 3, from which material parameters cy, C and n can be determined when employing the functional law extension (Eq. 2), must be narrowed between the initial instabilities 9 « 0.038 and strain at onset of necking 9 « 0.173. Material parameters, obtained by LSM are: cy = -73.517 MPa, C=744.68 MPa and «=0.1618. Note the negative value of cy which is unusual and in fact physically unacceptable. However, aiming to achieve the best fit along the monotonically increasing part of the yield curve the negative value of c y is obtained when excluding the beginning part of the data containing initial instabilities from the LSM procedure. Finally, the yield curve applicable for numerical simulation (Fig. 4) is assembled from the part containing initial instabilities, for 9<0.173 and following Eqs. 1, and functional law part (Eq. 2) which can be, in principle, extended beyond all strain limits. Inverse identification of the yield curve is based on a numerical simulation of the real tensile experiment. In fact it is achieved through a series of simulations in which, by appropriate tuning of the yield curve parameters that enter the numerical model, discrepancy between measured and computed tensile force is minimised. Since the yield curve (Fig. 3) before the onset of necking is considered reliable, only the post-necking part of the yield curve is sought. Thus, only the part of the force plot beyond the maximum measured force (Fig. 2) is intended for the identification purposes. The inverse methodology used for the yield curve characterization is explained in [10]; in the present paper only the result of the identification, 1.e. the yield curve, is given (Fig. 4). Comparing the two yield curves, derived respectively by Ludwig's law and inverse identification, two differences can be clearly seen in the post-necking part. First, slopes of the curves differ considerably. This difference is not always so obvious; namely, it depends on the tested material and its damage properties which are crucial for the material behaviour at large strains. The reason for different slopes is just in the consideration or no consideration of some 0.0 0.3 0.6 0.9 1.2 1.5 9 (/) Fig. 4. Yield curves by Ludwig's law and inverse identification damage mechanism. Since the background of all functional law yield curve extensions is to predict behaviour of the material after onset of necking just on the basis of the pre-necking data, it can not be expected that such a prediction is trustworthy. Hence, such a curve extension beyond the necking occurrence is a pure speculation. Contrary to this approach, the inverse identification includes available data of the post-necking behaviour into a characterization procedure. In those data the influence of damage mechanism is comprised. As a result, the yield curve is obtained in which material degradation is reflected. It is worth mentioning that a real damage mechanism is not recognized in this way, but its integral effect is nevertheless incorporated into the yield curve. Taking into consideration damage mechanisms with all their peculiarities (void nucleation and growth, stress triaxiality, etc.) would require complete constitutive modelling with new material laws, which are usually difficult to implement into commercial computer codes. Nowadays, an extensive research on various material damage phenomena and appropriate numerical implementation is in course [11] and [12], but it is beyond the subject of the present paper. The described approach with the yield curve tuning works fairly well for loading cases where material is loaded dominantly in tension, for example in deep drawing simulations. The second difference between the two considered yield curves is the existence of the termination point. While the functional law curve has no limit imposed on strain, the inversely Fig. 5. Tensile test set-up and modelling domain identified curve has a plain limit, which is based on the available post-necking data. Therefore no additional laws regarding material failure is needed. 2 YIELD CURVE VERIFICATION Both in the previous section identified yield curves (Fig. 4) will be used in computer simulations where material behaviour at large strains is of the prime concern, e.g. the necking formation at standard tensile test. Two simulations based on the same finite element model were performed, the only difference between the simulations being the adopted yield curve. Technical details regarding numerical simulation of such a test are given in [10]; here only explanation on the determination of a modelling domain is given. Namely, taking into account that necking of a tensile specimen usually develops in its central region, one can reduce the modelling domain at least to the part of the specimen between the extensometer clamps. Accordingly, appropriate boundary conditions must be set at this model boundary. Further, due to favourable prismatic shape of the reduced modelling domain three-fold symmetry can be introduced. Thus, only one eighth of the specimen between the extensometer clamps is modelled, as it is schematically shown in Fig. 5. Eight node brick elements, in finite element code ABAQUS [13] designated as C3D8I, are used for a domain discretization. Loading of the numerical model is achieved by applying a uniform displacement of all nodes on the extensometer clamp boundary plane (denoted with "A" in Fig. 5) in the longitudinal direction of the specimen. In this way stretching of the specimen in the same amount as measured in the real tensile test is assured. The calculated tensile force Fc, resulting from stretching of the numerical model, is obtained as a sum of nodal forces of the nodes laying in the opposite symmetry plane ("B" in Fig. 5). The resulting forces, obtained in the tensile test simulations by considering the considered two yield curves, are compared with the measured force Fm in Fig. 6. It can be seen that both calculated forces match perfectly the measured force till the necking occurrence at Fmax. This agreement confirms appropriate modelling of the experiment. In the continuation of the force plot, Fc obtained from the Ludwig's yield curve model significantly differs from Fm, while Fc obtained from the curve identified by the inverse method, matches Fm sufficiently well. The observed discrepancy is a direct consequence of considering different yield curves. In Figs. 7 and 8 the stress state and deformation in the necking region of the two models are shown, respectively. Fig. 7 represents Mises equivalent stress at the end of numerical simulation for the case of the inversely identified ALm (mm) Fig. 6. Measured and calculated force at tensile test Fig. 7. Mises equivalentMises in stresses region AL=28 mm - inversely identified yield curve yield curve. At the end of the simulation the model extension is AL=28 mm, which is 0.2 mm less than the extension measured in the real experiment. The final part of the yield curve used in the simulation has a steep descent (Fig. 4), which represents very fast softening, in fact tearing of the material. So, the stress value of 420 MPa, which is seen in Fig. 7, means almost total degradation (failure) of the material in the central part of the specimen. The shape (contour) of the neck region as well as thickness change of the numerical model are similar to those of the real specimen, which indicates proper modelling. Fig. 8 represents the deformation pattern at AL=28 mm (amount of the final stretch of the tensile specimen in the real experiment) for the case of Ludwig's yield curve. Excessive thinning of the necking region can be seen which is due to the chosen material constitutive model. Unrealistic large equivalent strain (9^5.4) is obtained solely because of monotonic hardening which is included in Ludwig's law by definition. In the shown example it is obvious that the response of the numerical model is wrong. In such a case the analyst must trace back through the iteration history of the simulation to find a believable stress-strain state. But using Ludwig's law no certain measure exists to find out, what is acceptable mechanical state after the onset of necking. Therefore, we can conclude that usage of pure functional extensional law, e.g. Ludwig's law, is anappropriate for simulations where large strains are expected. Fig. 8. Unrealistic deformation of necking region AL=28 mm - Ludwig's yield curve 3 DEEP DRAWING SIMULATION EXAMPLE In the previous example of a tensile test simulation the failure of correct modelling is clearly seen as an unrealistic mechanical state in the model. But, there exist cases where such inappropriateness is hardly to discover; an example is given in the sequel. In Fig. 9 a scheme of a classical deep drawing tool is shown. The tool consists of a fixed die, movable punch with stroke of H=260 mm and blankholder, that is pressed against the blank and die with a constant force of 6 MN (approximately 600 tons). Because of the two-fold symmetry of the drawn part, which is shaped as a square box of length 490 mm, width 387 mm, brim 140 mm and approximate drawn height of 260 mm with rounded corners, and considering also the orthotropic material, only V of the tool and blank is modelled. While the die and punch are in the numerical model taken as 0 50 100 150 200 250 300 Hpunch (mm) Fig. 10. Drawing force rigid bodies, blankholder is modelled as a deformable shell to account for a deformation of the drawing tools and press [14]. The blank is discretized with 3444 shell elements, in [13] designated as S4R. The material of the blank is the same as in the previous example of the tensile test simulations. Similarly, as in the tensile test example, two simulations will be conducted with the only difference between them being the yield curve used: the inversely identified vis-à-vis Ludwig's law yield curve. In Fig. 10 the drawing force executed by the punch during the forming process is plotted for both simulations. Because of the explicit formulation used in time (load) increment scheme both curves are significantly jagged; therefore, they need to be filtered to obtain readable data. One can observe coincidence of more than a half of the force plots, which is however expectable considering the almost equivalent beginning parts of their respective yield curves. In the continuation the force curves separate but not so Fig. 11. Mises equivalent stress - H=254 mm inversely identified yield curve much as in the tensile test example (Fig. 6). The reason is that the drawing force is a result of several factors: apart from the resistance exerted by the material of the blank with regard to the deformation imposed by the shape of the punch and die together, also friction between the blank and tools influence the stress distribution and shell thickness in critical cross-sections. Therefore the manifestation of the yield curve influence is smaller than in the tensile test, where no other factors except material strength govern the mechanical response of the tensile specimen. In Figs. 11 and 12 a comparison of Mises equivalent stress field for the mid-thickness layer of the shell at punch stroke of 254 mm is shown. Significantly different stress state can be seen. Due to large strains in case of the inversely identified yield curve, softening of the material initiates in the corner of the box. By further drawing tearing of the material is very likely to commence. This situation is shown in Fig. 11 where the white coloured area, denoted as "fracture", represents the domain of the box where tearing of the material is predicted. Since no finite element deleting criterion in case of fracture is employed in the present model, finite elements in the model still exist in the "fracture" area, but they have deteriorated stiffness properties. Presence of the material tearing can be perceived also in Fig. 10, where a sudden bend of the force curve appears in its final part. In Fig. 12, on the contrary, the softening phenomenon can not be seen since Ludwig's law has no softening. Stresses in deep drawing processes are to a great extent governed by stretching of the material, which is governed in turn by the tool shape and amount of the punch stroke. Since in both cases (Fig. 11 and 12) the punch stroke is equal and consequently the amount of drawing is similar, an Fig. 12. Equivalent Mises stress - H=254 mm Ludwig's yield curve average deformation is approximately the same, but stresses in Fig. 12 are, compared to Fig. 11, significantly lower. This is, of course, expectable due to the lower strength of the material, approximated by Ludwig's law. A comparison of calculated thicknesses of the drawn sheet is given in Figs. 13 and 14. Excessive thinning of the sheet can be seen in both cases. But, for the inversely identified yield curve case the thinning in Fig. 13 is accompanied with the softening in Fig. 11, thus the analyst can clearly see the potential danger. In the case of Ludwig's law yield curve no similar support from the stress state (Fig. 12) could be expected. Therefore the analyst can hardly forecast possible failure of the product without having additional information (experience from the field, similarities with an existing product, allowable thinning ratio, etc.). Remark: The steel material described in the present paper was not used in a production of the box, shown in the deep drawing example. This particular material and given deep drawing process parameters were selected only because the effects of two different yield curve characterization approaches are clearly seen. Under some other choices of material the crucial estimator which distinguishes between acceptable and inappropriate numerical modelling could be even harder to extract. 4 CONCLUSIONS In the paper a detail of a numerical simulation procedure which is important only at large strain state is discussed. Namely, at large strains, typically beyond the necking strain of the tensile test, often used approaches that rely on a functional law extension of the measured yield sth (mm) 1.50 1.25 1.00 0.75 0.50 0.25 0.00 curve, e.g. Ludwig's law, are improper. As an alternative, an inverse identification method of the yield curve characterization in the large strain region is proposed. Both approaches are examined in the tensile test simulation example, which serves excellently to show the deficiency of the functional law yield curve. In the tensile test simulation the response is completely governed by the yield curve and differences between the results of the two approaches are easy visible. The numerical model built by considering Ludwig's law yield curve fails to reproduce the tensile test measured response. In the deep drawing simulation example differences between the two approaches are still remarkable. While the inverse identification approach forecasts tearing of the drawn part, the Ludwig's law approach predicts integrity of the part. The only sign of possible fracture is very thin shell in the corner region of the box, but without knowing what is allowable thinning, no sure prediction can be given. Thus, the deep drawing example shows the real problem of inappropriate yield curve modelling. Namely, one can conduct a simulation of the high straining process with inappropriate yield curve and, if no support in the form of additional information is available, wrong prediction of the forming process, based on the misleading results of the simulation, is inevitable. This is in fact the main point which we want to highlight in the present paper. In the proposed inverse identification of the yield curve a classical metal plasticity model which employs Hill's yield surface is used. No additional constitutive law, where damage of the material is addressed in a more appropriate way, as for example in the Gurson-Twergaard- 0.23 sth (mm) •1.50 •1.25 •1.00 ■0.75 ■0.50 ■0.25 ■0.00 Fig. 13. Shell thickness - H=254 mm inversely identified yield curve Fig. 14. Shell thickness -H=254 mm Ludwig's yield curve Needleman (GTN) model, is used. Since in the proposed inverse identification damage of the material is enveloped by the yield curve which is characterized on the basis of the tensile test, also usefulness of the identified yield curve is limited to the material stretching cases, deep drawing processes being surely a part of them. ACKNOWLEDGEMENTS The authors gratefully acknowledge the staff at the laboratory for measurement of the Kovinoplastika Lož company for their excellent performance of all the measurements, referred to in the paper. 5 REFERENCES [1] European Committee for Standardization (1990) European Standard EN 10 002-1: Metallic materials - Tensile testing, Part 1: Method of test at ambient temperature. [2] Gerlach J., Keßler L. Material Parameters for the FEM-Simulation - Impact of Advanced Material Models on Forming and Part Properties. In: Kuzman K. et all. (Eds.): Proceedings of the IDDRG 2003 Conference, Bled, p. 97-107, 2003. [3] Gusek C.O., Bleck W., Dahl W. Modelling of sheet metal testing, Computational Materials Science, 1996 (7), p. 173-180. [4] Hambli R., Badie-Levet D. FEM Characterization of damage Law Parameters by Inverse Technique. In: Shirvani B. et all. (Eds.): Sheet Metal 2000, Proceedings of the 8th International Conference, p. 511-518, 2000. [5] AutoForm, AutoForm Engineering GmbH, 2008. [6] PAM-STAMP 2G, ESI Group, 2008. [7] Brigman P.W. Studies in Large Plastic Flow and Fracture, McGraw-Hill, 1952, New York. [8] Ling Y. Uniaxial True Stress-Strain after Necking, AMP Journal of Technology, 1996 (5), p. 37-48. [9] Sol H., Oomens C.W.J. Proceedings of the Euromech Colloquium held in Kerkrade, Kluwer Academic Publishers, Dordrecht, 1997. [10] Koc P., Štok B. Computer-aided identification of the yield curve of a sheet metal after onset of necking, Computational Materials Science 2004 (31), p. 155-168. [11] Halilovič, M., Vrh, M., Štok, B. Prediction of elastic recovery of a formed steel sheet considering stiffness degradation. Meccanica, (in press). [12] Vrh, M., Halilovič, M., Štok, B. Impact of Young's modulus degradation on springback calculation in steel sheet drawing, Strojniški vestnik - Journal of Mechanical Engineering, 2008 (54), p. 288-296. [13] ABAQUS Version 6.7, ABAQUS, Inc. and Dassault Systèmes, 2008. [14] Mole, N., Poje, J., Štok, B. Effect of selective lubrication and variable blank holding on deep drawing conditions. In: Juster, N. P., Rosochowski, A. (Eds.): The 9th International ESAFORM Conference on Material Forming, Glasgow, p. 255-258, 2006.