Prediction of Wrinkling of a Beverage Can Subjected to the Redrawing Process by J2 Deformation Theory

: Wrinkling of beverage cans is one of the problems faced by can manufacturers and aluminum suppliers. The bottom of an aluminum can is wrinkled by compression during the forming process. In this study, to predict the occurrence of wrinkles during the redrawing process of AA3104 ( t = 0.265 mm), which is the material used to fabricate aluminum cans, the classical plasticity J2 deformation theory (J 2 D) and ﬂow theory (J 2 F) were considered. J 2 F considers only the deformation perpendicular to the yield locus, whereas J 2 D considers the deformation perpendicular to the yield locus and that tangential to the yield locus. Wrinkles are predicted using ﬁnite element (FE) analyses based on J 2 D and J 2 F, and the results are compared. J 2 F could not predict the number and amplitude of wrinkles. By contrast, the wrinkles predicted using J 2 D exhibited good agreement with sample data obtained for a real can. To ﬁnd the di ﬀ erence between the results obtained using J 2 F and J 2 D, evolutions of stress path in a wrinkled element are compared. It was conﬁrmed that compressive stress is more dominant in the J 2 D case than in the J 2 F case. Moreover, the measured e ﬀ ective strain of the element is small, under 0.04. In conclusion, J 2 D is more suitable for predicting the wrinkling behavior of aluminum cans than J 2 F. In addition, ANOVA and ANOM analysis are performed to evaluate the inﬂuence of the design parameters, namely friction coe ﬃ cient, thickness, and outer proﬁle angle, and the parameters are optimized to reduce wrinkles by combining the Taguchi method with FE simulation based on the J 2 D theory.


Introduction
Wrinkling problem is encountered in the drawing of automotive parts and beverage cans for over 40 years. A known cause of wrinkles is compressive stress due to plastic behavior during the formation of thin sheets. Indeed, the aluminum manufacturers, can makers manufacturing techniques and cost control efforts have resulted in thinner sheet thicknesses, thereby causing wrinkle problems. To avoid wrinkling and obtain the desired shape, the effect of the blank holding force, clearance between punch and lower die, tool design, and lubrication have been studied. When the blank holding force is increased to reduce wrinkling, problems, such as thinning and fracture of the sheet, should be considered. In addition, wrinkling is affected by many factors, such as material properties and thickness, and tool geometry, such as outer profile angle, coefficient of friction, and lubrication. Especially, these factors are more sensitive in the case of beverage cans because the blank is very thin, ranging from 0.2 to 0.4 mm. Moreover, the process of aluminum beverage can is consists of drawing and ironing (D&I). The wrinkles in the beverage can are known to occur during the redrawing process. Figure 1 shows a schematic diagram of the process to fabricate beverage cans. The occurrence of wrinkles on the bottom of the beverage can is shown in Figure 2. In the same figure, many small wrinkles can be observed. To overcome such wrinkling, optimization approaches have been investigated in previous studies [1][2][3]. Moreover, process designers and engineers are expanding the application of finite element (FE) simulations [4,5]. FE simulations have been used to predict wrinkling for more than 30 years now [6][7][8][9][10][11]. Cao et al. [12][13][14][15] employed FE simulation to predict wrinkling behavior in sheet metal forming. Makinouchi [16] conducted FE simulations to study the effect of blank holding on wrinkling in automotive parts. Nam and Han [17] studied the anisotropic effects of low-carbon steel sheets by comparing their thickness distributions after D&I. Rekas [18] studied the effect of the thickness of AA3XXX grade sheet on wrinkling by performing an FE simulation. The author identified the optimum sheet thickness for obtaining a good product given that the proper tool geometry, process parameters, and clearance are selected. Kawka et al. [19] predicted the wrinkles formed after deep drawing. They considered plastic behaviors by using the hardening function and the anisotropic yield function to describe material properties accurately. Kim and Son [20] approached wrinkling in terms of various material parameters, such as yield stress, hardening coefficient, and perpendicular anisotropy. They reported that wrinkles are formed rapidly as the hardening exponent n decreases, initial yield stress increases, and thickness decreases. Moreover, they investigated the wrinkle limits numerically. Recently, the wrinkle prediction problem was presented as a benchmark at NUMISHEET 2014. Kim et al. [21] predicted the onset and growth of wrinkles by adopting a bifurcation algorithm to implement cylindrical cup deep drawing. Most researchers have pointed out that tool path, die, blank holder, heating condition, friction, and appropriate material properties must be considered simultaneously to obtain accurate FE simulation results.    Figure 1 shows a schematic diagram of the process to fabricate beverage cans. The occurrence of wrinkles on the bottom of the beverage can is shown in Figure 2. In the same figure, many small wrinkles can be observed. To overcome such wrinkling, optimization approaches have been investigated in previous studies [1][2][3]. Moreover, process designers and engineers are expanding the application of finite element (FE) simulations [4,5]. FE simulations have been used to predict wrinkling for more than 30 years now [6][7][8][9][10][11]. Cao et al. [12][13][14][15] employed FE simulation to predict wrinkling behavior in sheet metal forming. Makinouchi [16] conducted FE simulations to study the effect of blank holding on wrinkling in automotive parts. Nam and Han [17] studied the anisotropic effects of low-carbon steel sheets by comparing their thickness distributions after D&I. Rekas [18] studied the effect of the thickness of AA3XXX grade sheet on wrinkling by performing an FE simulation. The author identified the optimum sheet thickness for obtaining a good product given that the proper tool geometry, process parameters, and clearance are selected. Kawka et al. [19] predicted the wrinkles formed after deep drawing. They considered plastic behaviors by using the hardening function and the anisotropic yield function to describe material properties accurately. Kim and Son [20] approached wrinkling in terms of various material parameters, such as yield stress, hardening coefficient, and perpendicular anisotropy. They reported that wrinkles are formed rapidly as the hardening exponent n decreases, initial yield stress increases, and thickness decreases. Moreover, they investigated the wrinkle limits numerically. Recently, the wrinkle prediction problem was presented as a benchmark at NUMISHEET 2014. Kim et al. [21] predicted the onset and growth of wrinkles by adopting a bifurcation algorithm to implement cylindrical cup deep drawing. Most researchers have pointed out that tool path, die, blank holder, heating condition, friction, and appropriate material properties must be considered simultaneously to obtain accurate FE simulation results.   Under the classical plasticity theory, two theories have been developed to calculate strain, namely flow theory and deformation theory. The two theories differ in the manners in which they compute strain rate. Flow theory has been used widely because it consists of simple governing equation terms and it is easy to implement in FE simulation programs. Hutchinson et al. [22][23][24][25][26] studied and compared deformation theory and flow theory by using the von Mises yield function based on an invariant J2. Cao et al. [27] showed that a wide range of material characteristics can be applied by implementing FE simulation based on the kinematic-hardening model and the J2 flow theory (J 2 F) or J2 deformation theory (J 2 D).
Neale and Hassani [28] reported the influences of J 2 D and J 2 F on the wrinkling limit diagram of an isotropic material. They found that the deformation theory of plasticity predicted wrinkling at lower stress levels than the flow theory, as expected, and the deformation theory generated more conservative predictions of the onset of wrinkling.
Dick and Yoon [29] studied the effects of the hardening function and anisotropic yield function from the viewpoint of considering the properties of anisotropic material behavior. Moreover, they performed an FE simulation by using the flow theory in conjunction with Barlat's Yld2000-2d anisotropic yield function [30], which provides a good description of anisotropic material behavior. To evaluate the hardening effect, they considered four different hardening curves of AA 5042 and compared wrinkling behaviors, punch forces, and punch displacements. The comparison results indicated that a lower yield stress significantly reduced the magnitude of wrinkling. Although the predicted punch force was higher than the experimental value, the results obtained using the four different hardening curves were almost identical. Moreover, the results of the simulation performed using the Yld2000-2d anisotropic yield function agreed well with the experimental results. The magnitude of the wrinkles was larger than 2 mm and the number of wrinkles was small, that is 13 wrinkles over 360 • .
In this study, we considered beverage cans because more than 40 wrinkles of small magnitude are created in the fabrication of this product. We performed FE simulations based on the classical plasticity theories of J 2 D and J 2 F to predict the occurrence of wrinkles. Moreover, we considered only the redrawing process because wrinkles are formed in this process. According to this consideration, the blank was assumed to have a homogeneous thickness of 0.265 mm because this thickness is obtained after deep drawing in a real process. To perform the FE simulation, commercial software ABAQUS/STANDARD was used for wrinkle prediction. In addition, Taguchi's orthogonal array simulations were combined with the FE simulations to quantify and optimize the main variables affecting the forming quality of a can. Three design variables, namely friction coefficient, blank thickness, and outer profile angle, were optimized. The simulation results obtained using the optimized variables indicated reduced wrinkle formation.

J2 Deformation Theory
Storen and Rice [31] proposed the J2 deformation theory based on the necking phenomenon, which is a local plastic deformation caused by the instability of plane stresses. This theory was derived by assuming a local bifurcation by using the difference in the velocity rates of deformation between the interior and exterior sides of a thin sheet. In J 2 D, the hypothesis is initial isotropic in plastic behavior and plastic deformation covers only changes in shape, not changes in volume. The second invariant J2 was employed as the isotropic yield function, commonly called the von Mises yield function. The yield function can be expressed as follows: where σ 1 , σ 2 , and σ 3 denote the principal stresses. In the case of plane stress, the above function can be simplified as follows: (2) J 2 D can be expressed in terms of equivalent stress (σ) and equivalent plastic strain (ε p ). The equivalent stress and equivalent plastic strain can be written as follows: here σ ij and ε p ij denote the deviatoric stress and plastic strain tensor, respectively. By using Hencky's equation, one can express the plastic strain (ε p ij ) as follows: where ∅ and ε p denote a proportional constant and the plastic equivalent strain, respectively. Moreover, the plastic equivalent strain can be rewritten in terms of total strain (dε) and elastic strain (dε e ). E s = σ ε and E = σ ε e denote the secant modulus and modulus of elasticity obtained from the stress-strain curve, respectively.
Using the above equation, one can express the derivative term ( . ∅) of the proportional constant (∅) as follows: where E t denotes the tangential modulus of elasticity, and . σ represents the velocity for equivalent stress. The velocity of plastic strain (ε p ij ) in Equation (4) can be written as follows: This equation can be expressed simply by substituting A and B as follows: (8) is the mathematical expression of the deformation theory.
In the above equation, the first and second terms on the right-hand side of Equation (8) describe the plastic deformation in the vertical direction and the tangential direction from the yield locus, respectively. In contrast to J 2 F, which assumes that plastic deformation occurs only perpendicular to the yield locus, J 2 D permits deformation in both the vertical and tangential directions. This can be confirmed from Equation (8).

J2 Flow Theory
J 2 F, too, uses an isotropic yield function based on the second invariant J2. In Equation (8), when E s = E, A = 0 and B = 9 E . Therefore, Equation (8) becomes the equation for J 2 F and it can be expressed as follows: From the above equation, it can be confirmed that J 2 F considers plastic deformation only in the perpendicular direction.

Uniaxial Tensile Test
To confirm the flow stress-strain relationship and the anisotropic behavior of an AA3104 sheet, uniaxial tensile tests were performed at θ = 0 • , 45 • , and 90 • relative to the rolling direction (RD). Specimens were prepared using a laser machine according to the ASTM-E8 standard. The tensile tests were carried out at a strain rate of approximately 1 mm/min. The initial yield stress and r-value were measured in the same test. Figure 3 shows the true stress-strain curves at θ = 0 • , 45 • , and 90 • relative to the RD. In the figure, the results of the tensile test depict isotropic material behavior, but the r-values (0 • : 0.41, 45 • : 0.99, and 90 • : 1.02) indicate anisotropic material behavior. The mechanical properties obtained from the uniaxial tension test are summarized in Table 1, including modulus of elasticity (E), Poisson's ratio (υ), ultimate tensile strength (UTS), initial yield stress (YS), and r-values.

Uniaxial Tensile Test
To confirm the flow stress-strain relationship and the anisotropic behavior of an AA3104 sheet, uniaxial tensile tests were performed at θ = 0°, 45°, and 90° relative to the rolling direction (RD). Specimens were prepared using a laser machine according to the ASTM-E8 standard. The tensile tests were carried out at a strain rate of approximately 1 mm/min. The initial yield stress and r-value were measured in the same test. Figure 3 shows the true stress-strain curves at θ = 0°, 45°, and 90° relative to the RD. In the figure, the results of the tensile test depict isotropic material behavior, but the r-values (0°: 0.41, 45°: 0.99, and 90°: 1.02) indicate anisotropic material behavior. The mechanical properties obtained from the uniaxial tension test are summarized in Table 1, including modulus of elasticity (E), Poisson's ratio (υ), ultimate tensile strength (UTS), initial yield stress (YS), and r-values.

Geometrical Model
Can manufacturing involves multiple deep drawing processes, including redrawing, three-step ironing, doming, trimming, and flanging. However, because wrinkles occur during redrawing, only the redrawing process was considered for wrinkle prediction. A schematic drawing of the redrawing process and the corresponding geometry for FE simulation along with a mesh of the blank are shown in Figure 4a,b, respectively. The FE simulation model consisted of a redraw die, blank holder, punch, cup locator, and blank with 1/4 symmetry. The tool and can dimensions were set by referencing the size of a commercially manufactured 12-oz (355 mL) D&I can [18,19]. In addition, the external profile angle was modeled as α = 37°. The redraw die and punch were set as a rigid body for analysis. The four-node shell element (S4R) and five integration points in the thickness direction were applied and the smallest element size of the blank was set to 1 mm × 0.9 mm. A blank holding force of approximately 11 kN was applied. Moreover, the friction coefficient between all parts (the blank holder and the plate, the punch and the plate, the cup locator and the blank) was applied in the form of penalty friction coefficient of 0.05 since the range of that is from 0.03 to 0.07 for a real can making process [17,29,32]. The punch stroke was set to be 13.5 mm from contact to the final step.

Geometrical Model
Can manufacturing involves multiple deep drawing processes, including redrawing, three-step ironing, doming, trimming, and flanging. However, because wrinkles occur during redrawing, only the redrawing process was considered for wrinkle prediction. A schematic drawing of the redrawing process and the corresponding geometry for FE simulation along with a mesh of the blank are shown in Figure 4a,b, respectively. The FE simulation model consisted of a redraw die, blank holder, punch, cup locator, and blank with 1/4 symmetry. The tool and can dimensions were set by referencing the size of a commercially manufactured 12-oz (355 mL) D&I can [18,19]. In addition, the external profile angle was modeled as α = 37 • . The redraw die and punch were set as a rigid body for analysis. The four-node shell element (S4R) and five integration points in the thickness direction were applied and the smallest element size of the blank was set to 1 mm × 0.9 mm. A blank holding force of approximately 11 kN was applied. Moreover, the friction coefficient between all parts (the blank holder and the plate, the punch and the plate, the cup locator and the blank) was applied in the form of penalty friction coefficient of 0.05 since the range of that is from 0.03 to 0.07 for a real can making process [17,29,32]. The punch stroke was set to be 13.5 mm from contact to the final step.

Material Model
FE simulations based on J2F and J2D were performed using commercial software ABAQUS/STANDARD. To describe the material behavior during the redrawing process, the Ramberg-Osgood [33] plasticity model was adopted to represent the flow stress of AA3104 material in the FE simulation. A one-dimensional model of this equation can be expressed as follows: where E, ε, σ, K and n are elastic modulus, strain, stress, material constant, and hardening exponent, respectively. By introducing the yield strength of the material, , and defining a new parameter, α, related to K as = ( / ) , the Equation (10) can be rewritten as follows: here, E = 69.73 GPa, = 281 MPa, α = 0.45, and n = 22.59 were applied to explain the flow curve of AA3104. The curve obtained using Equation (11) and experimental data (0°) are shown in Figure 5. The result of the fitting curve shows that the Ramberg-Osgood model can describe the experiment well.

Material Model
FE simulations based on J 2 F and J 2 D were performed using commercial software ABAQUS/STANDARD. To describe the material behavior during the redrawing process, the Ramberg-Osgood [33] plasticity model was adopted to represent the flow stress of AA3104 material in the FE simulation. A one-dimensional model of this equation can be expressed as follows: where E, ε, σ, K and n are elastic modulus, strain, stress, material constant, and hardening exponent, respectively. By introducing the yield strength of the material, σ 0 , and defining a new parameter, α, related to K as α = K(σ 0 /E) n−1 , the Equation (10) can be rewritten as follows: here, E = 69.73 GPa, σ 0 = 281 MPa, α = 0.45, and n = 22.59 were applied to explain the flow curve of AA3104. The curve obtained using Equation (11) and experimental data (0 • ) are shown in Figure 5. The result of the fitting curve shows that the Ramberg-Osgood model can describe the experiment well.

= +
here, E = 69.73 GPa, = 281 MPa, α = 0.45, and n = 22.59 were applied to explain the flow curve of AA3104. The curve obtained using Equation (11) and experimental data (0°) are shown in Figure 5. The result of the fitting curve shows that the Ramberg-Osgood model can describe the experiment well.

Wrinkling Prediction
The strain distribution at the final stage was analyzed with both J 2 F and J 2 D, as shown in Figure 6a,b, respectively. In this figure, the strain distribution is uniform along the radial direction because it is based on the second invariant J2. In addition, the wrinkled sample used in this study was obtained from a can manufacturer and the wrinkles were measured at z = 3.5 mm from the bottom of the can because the maximum wrinkling was observed at this point. To compare the wrinkles indicated by each theory and the actual wrinkles on the sample can, the radial distance was measured at the same position (z = 3.5 mm) on the sample can (Figure 6c). A schematic diagram of the wrinkle measurement process is illustrated in Figure 6d. ∆h denotes the highest wrinkle amplitude and can be obtained from ∆h = r max − r min . Only one sample was used since the characteristics, which are amplitudes and the number of wrinkles, are similar in most samples. During experimental processes, only one ∆h was determined for one sample and each FE model also provided just one value of ∆h. In this study, the value of ∆h was determined directly from the measurement without any deviation because the wrinkles were measured for only one sample.

Wrinkling Prediction
The strain distribution at the final stage was analyzed with both J2F and J2D, as shown in Figure  6a,b, respectively. In this figure, the strain distribution is uniform along the radial direction because it is based on the second invariant J2. In addition, the wrinkled sample used in this study was obtained from a can manufacturer and the wrinkles were measured at z = 3.5 mm from the bottom of the can because the maximum wrinkling was observed at this point. To compare the wrinkles indicated by each theory and the actual wrinkles on the sample can, the radial distance was measured at the same position (z = 3.5 mm) on the sample can (Figure 6c). A schematic diagram of the wrinkle measurement process is illustrated in Figure 6d. ∆ℎ denotes the highest wrinkle amplitude and can be obtained from ∆ℎ = − . Only one sample was used since the characteristics, which are amplitudes and the number of wrinkles, are similar in most samples. During experimental processes, only one ∆h was determined for one sample and each FE model also provided just one value of ∆ℎ. In this study, the value of ∆ℎ was determined directly from the measurement without any deviation because the wrinkles were measured for only one sample.
Eleven wrinkles were observed on the sample obtained from the manufacturer. However, four wrinkles were predicted in the simulation based on J2F, while 10 wrinkles were predicted in the simulation based on J2D. The ∆ℎ, which means the maximum value of amplitudes measured on the sample can, was approximately 0.108 mm and the simulations based on J2F and J2D yielded ∆ℎ predictions of approximately 0.009 mm and 0.077 mm, respectively. To identify the evolution of the stress path at the point of wrinkle occurrence, the stress path obtained from an element in which wrinkles occurred and the initial von Mises yield locus are shown in Figure 7. The circle represents the initial yield locus, and the solid and dotted lines show the stress paths obtained in the simulations based on J2D and J2F, respectively. In addition, the pure shear line is shown to confirm the dominance of the compressive stress region, where the left side of the pure Eleven wrinkles were observed on the sample obtained from the manufacturer. However, four wrinkles were predicted in the simulation based on J 2 F, while 10 wrinkles were predicted in the simulation based on J 2 D. The ∆h, which means the maximum value of amplitudes measured on the sample can, was approximately 0.108 mm and the simulations based on J 2 F and J 2 D yielded ∆h predictions of approximately 0.009 mm and 0.077 mm, respectively.
To identify the evolution of the stress path at the point of wrinkle occurrence, the stress path obtained from an element in which wrinkles occurred and the initial von Mises yield locus are shown in Figure 7. The circle represents the initial yield locus, and the solid and dotted lines show the stress paths obtained in the simulations based on J 2 D and J 2 F, respectively. In addition, the pure shear line is shown to confirm the dominance of the compressive stress region, where the left side of the pure shear line is predominant. Both stress paths are similar in the elastic region. However, in the plastic region, the stress path measured from J 2 F is shifted to the upper of the pure shear line (second quadrant), while that obtained from J 2 D is shifted toward the left of the pure shear line (second quadrant). This means that the stress path of J 2 D is in the region dominated by compressive stress compared to that of J 2 F. Because wrinkle formation is dominated by compressive stress, the prediction obtained using J 2 F is superior to that obtained using J 2 D. This result is consistent with the results of previous studies, which indicate that wrinkle formation is governed by compressive stress. The magnitudes of stresses and strains for both cases in the final stroke (13.  . This means that the stress path of J2D is in the region dominated by compressive stress compared to that of J2F. Because wrinkle formation is dominated by compressive stress, the prediction obtained using J2F is superior to that obtained using J2D. This result is consistent with the results of previous studies, which indicate that wrinkle formation is governed by compressive stress. The magnitudes of stresses and strains for both cases in the final stroke (13.5 mm) were measured: = −126 MPa, = 199.8 MPa and ε = −0.027, ε = 0.026 for J2F and = −248.0 MPa, = 83.8 MPa and ε = −0.034, ε = 0.024 for J2D. Moreover, the stress and strain values indicate that J2D is in the area of plastic behavior, where compression is predominant.
To verify the wrinkling behavior, the stresses at a cylindrical coordinate, which uses the radial, circular, and axial directions, were analyzed for J2F and J2D. Figures 8 and 9 show the predictions of stress distributions on the outer surface of the deformed specimen based on the J2F and J2D theories, respectively. Both figures indicate a negative value of and positive value of . In other words, the material was stretched in the radial direction and compressed in the circular direction, which can be seen in the figures. The stress distributions were measured at punch stroke = 0, 5, wrinkle initiation stroke, and final stroke. In both figures, the punch strokes of 0 and 5 mm did not lead to wrinkle formation, while the wrinkle initiation stroke (10.97 mm for J2F and 11.36 mm for J2D) and the final stroke led to stress distributions with wrinkle formation. However, J2F predicted the formation of a small number of wrinkles, whereas J2D predicted the formation of a large number of wrinkles. Moreover, the stretched stress and compressive stress at the initiation stroke were more dominant according to J2D than J2F, but those of them in the final stroke show converse results. The radial stress in the final stroke was more dominant according to J2F than J2D. The stresses at the initiation stroke and final stroke were as follows:  [28]. Moreover, the difference between global and cylindrical coordinate system stresses comes from coordinate systems.  To verify the wrinkling behavior, the stresses at a cylindrical coordinate, which uses the radial, circular, and axial directions, were analyzed for J 2 F and J 2 D. Figures 8 and 9 show the predictions of stress distributions on the outer surface of the deformed specimen based on the J 2 F and J 2 D theories, respectively. Both figures indicate a negative value of σ θ and positive value of σ r . In other words, the material was stretched in the radial direction and compressed in the circular direction, which can be seen in the figures. The stress distributions were measured at punch stroke = 0, 5, wrinkle initiation stroke, and final stroke. In both figures, the punch strokes of 0 and 5 mm did not lead to wrinkle formation, while the wrinkle initiation stroke (10.97 mm for J 2 F and 11.36 mm for J 2 D) and the final stroke led to stress distributions with wrinkle formation. However, J 2 F predicted the formation of a small number of wrinkles, whereas J 2 D predicted the formation of a large number of wrinkles. Moreover, the stretched stress and compressive stress at the initiation stroke were more dominant according to J 2 D than J 2 F, but those of them in the final stroke show converse results. The radial stress in the final stroke was more dominant according to J 2 F than J 2 D. The stresses at the initiation stroke and final stroke were as follows: J 2 F (σ θ = −185.13 MPa, σ r = 182.88 MPa for wrinkle initiation, and σ θ = −217.84 MPa, σ r = 237.11 MPa for final stroke); J 2 D (σ θ = −222.97 MPa, σ r = 223.56 MPa for wrinkle initiation and σ θ = −216.27 MPa, σ r = 217.67 MPa for final stroke), respectively. These results are consistent with the results of Neale and Tugcu [28]. Moreover, the difference between global and cylindrical coordinate system stresses comes from coordinate systems. To identify the cause of wrinkle occurrence, we observed the contact region between the tool and the blank. Figure 10a shows the redrawing tool and blank of half of the model. In the figure, the tool is not contacting the blank in the wrinkled area during redrawing. It is easy to apply compressive stress to the blank. Because of this phenomenon, the blank is wrinkled easily. Moreover, the stress mode of a wrinkled element is illustrated in Figure 10b to depict the wrinkling mechanism. On the To identify the cause of wrinkle occurrence, we observed the contact region between the tool and the blank. Figure 10a shows the redrawing tool and blank of half of the model. In the figure, the tool is not contacting the blank in the wrinkled area during redrawing. It is easy to apply compressive stress to the blank. Because of this phenomenon, the blank is wrinkled easily. Moreover, the stress To identify the cause of wrinkle occurrence, we observed the contact region between the tool and the blank. Figure 10a shows the redrawing tool and blank of half of the model. In the figure, the tool is not contacting the blank in the wrinkled area during redrawing. It is easy to apply compressive stress to the blank. Because of this phenomenon, the blank is wrinkled easily. Moreover, the stress mode of a wrinkled element is illustrated in Figure 10b to depict the wrinkling mechanism. On the cylindrical coordinates, the element is subjected to tensile stress along the radial direction ®, which is indicative of the drawing mode. In contrast, it is subjected to compressive stress along the angle direction (θ). Wrinkles are easily generated in this mode. Because wrinkles are formed according to this mechanism, the outer profile angle, which affects to condition of contact, is an important design parameter. Moreover, to confirm the strain value at the element in which wrinkling occurred, the effective strains according to the punch stroke were determined using the J 2 D and J 2 F theories and compared with each other (Figure 11). Both strains were less than 0.04. This result indicated that a small deformation occurred in the wrinkled area. Based on the above observations, J 2 D is more suitable for wrinkle prediction than J 2 F.

Wrinkle Optimization
In the light of the above observations, Taguchi's orthogonal arrays were employed in conjunction with FE simulation based on J2D to optimize the redrawing process conditions to avoid wrinkle formation. In a Taguchi orthogonal array experiment, the objective function to be optimized must be defined. Moreover, the signal-to-noise (S/N) ratio, a higher value of which indicates a superior result, was used. The S/N ratio is defined as follows: where ηi and ∆hi are the observed value (in dB) and the highest amplitude at z = 3.5 mm from the bottom of the can in the i-th simulation, respectively. The highest value of the S/N ratio corresponds to the best condition.

Wrinkle Optimization
In the light of the above observations, Taguchi's orthogonal arrays were employed in conjunction with FE simulation based on J2D to optimize the redrawing process conditions to avoid wrinkle formation. In a Taguchi orthogonal array experiment, the objective function to be optimized must be defined. Moreover, the signal-to-noise (S/N) ratio, a higher value of which indicates a superior result, was used. The S/N ratio is defined as follows: where ηi and ∆hi are the observed value (in dB) and the highest amplitude at z = 3.5 mm from the bottom of the can in the i-th simulation, respectively. The highest value of the S/N ratio corresponds to the best condition. Moreover, FE simulation, according to the matrix and the evaluation of the S/N ratio with

Wrinkle Optimization
In the light of the above observations, Taguchi's orthogonal arrays were employed in conjunction with FE simulation based on J 2 D to optimize the redrawing process conditions to avoid wrinkle formation. In a Taguchi orthogonal array experiment, the objective function to be optimized must be defined. Moreover, the signal-to-noise (S/N) ratio, a higher value of which indicates a superior result, was used. The S/N ratio is defined as follows: where η i and ∆h i are the observed value (in dB) and the highest amplitude at z = 3.5 mm from the bottom of the can in the i-th simulation, respectively. The highest value of the S/N ratio corresponds to the best condition. Moreover, FE simulation, according to the matrix and the evaluation of the S/N ratio with analysis of the mean (ANOM), which determine the priority of influences and the analysis of the variance (ANOVA) for determining quality characteristics, afforded a better understanding of the individual effects of each factor. The factors and their levels in the Taguchi method are summarized in Table 2. Three factors, namely A-friction coefficient, B-can thickness and C-outer profile angle, and three levels, as well as the values of the three factors and each level, are represented. Based on Table 3, the results of a matrix simulation involving the L 9 orthogonal array selected from Taguchi's standard orthogonal array table, objective value (∆h), and S/N ratio are summarized in Table 4. The numbers in each column are related to the level number of each factor. ∆h is the value measured along the same lines as in Figure 6a. In this study, we investigated only the individual effect of each factor on the quality characteristics without considering the interactions among the factors. The FE simulation of the redrawing process was then performed according to the L 9 orthogonal array by using J 2 D. To evaluate the influence of each factor, the sum of the squares (SS) method was used. SS due to factor j (j = A, B, C) was calculated as follows: where m is the overall mean of the η i values in the simulations, defined as m = 1/9 η i = 40.56 (dB), where m ji is the average of η related to the level i (i = 1, 2, 3) of factor j given by m ji = 1/3 (η j ) i . ANOM and ANOVA were calculated from the matrix simulations along the L 9 orthogonal array. Moreover, the result of ANOM for the redrawing simulation is presented in Figure 12, which shows the effects of individual factors on the variation of the S/N ratio. The outer profile angle significantly affected wrinkle occurrence during redrawing and was a major control factor for optimizing the redrawing condition. Moreover, the S/N ratio increased when the friction coefficient increased, sheet thickness increased, and outer profile angle decreased. The results of ANOM and ANOVA are summarized in Table 4. The factor levels corresponding to the highest S/N ratio were selected to optimize the condition.
The optimum values of the factors and their levels are given in Table 4. The optimal values of each factor and level are indicated by an asterisks (*). From the table, the optimal condition for minimizing wrinkles is A 3 B 3 C 1 , which indicates that the S/N ratio of the friction coefficient µ is optimized in level 3 (µ = 0.07), blank thickness is optimized in level 3 (t = 0.28 mm), and outer profile angle is optimized in level 1 (α = 34). Moreover, according to the ANOVA result, the outer profile angle had the strongest influence on wrinkle occurrence.   To verify the optimization results, an FE simulation was carried out under the optimal condition A 3 B 3 C 1 . Figure 13 shows the results of this simulation and the wrinkle amplitudes along the circumferential direction. The results indicated that the number and amplitude of wrinkles were minimized. The measured amplitude of the wrinkles was 0.0001 mm.

Conclusions
FE simulation is widely used to predict wrinkling behavior in redrawing for manufacturing beverage cans. In this study, an FE simulation of the redrawing process was performed with two classical plasticity theories, namely J2F and J2D, which are based on the second invariant J2. The

Conclusions
FE simulation is widely used to predict wrinkling behavior in redrawing for manufacturing beverage cans. In this study, an FE simulation of the redrawing process was performed with two classical plasticity theories, namely J 2 F and J 2 D, which are based on the second invariant J2. The number and amplitude of wrinkles in a beverage can that was redrawn using AA3104 material (t = 0.265 mm) were predicted and compared to verify the effects of each theory. The number and amplitude of wrinkles were measured at the same height on the sample in simulations conducted using the J 2 F and J 2 D theories. In conclusion, we demonstrated the feasibility of using the J 2 D theory to predict the formation of many small-amplitude wrinkles in the FE simulation. The results can be summarized as follows: 1.
The equivalent plastic strain in the redrawing process for manufacturing beverage cans was as small as <0.04.

2.
J 2 F underestimated the amplitude and number of wrinkles.

3.
The amplitude and the number of wrinkles predicted using J 2 D was shown to be in agreement with the data measured from the sample.

4.
The stress paths obtained using J 2 F and J 2 D were compared to confirm the difference between them and it was confirmed that the stress path obtained using J 2 D was in the region dominated by compression.

5.
To optimize the design parameters, the Taguchi method and ANOVA and ANOM analyses were employed along with FE simulations based on J 2 D. 6.
The optimal condition for minimizing wrinkle formation was A 3 B 3 C 1 , which means that the S/N ratio of the friction coefficient µ was optimized in level 3 (µ = 0.07), the blank thickness was optimized in level 3 (t = 0.28), and the outer profile angle was optimized in level 1 (α = 34), respectively. 7.
A simulation with the optimal condition shows the result without wrinkle.
Based on the above results, the validity of the J 2 D theory was demonstrated. When predicting many small-amplitude wrinkles, selection of the plastic strain calculation theory is more important than selection of the yield function. In addition, the optimization was performed by combining FE simulation based on J 2 D and the Taguchi method to indicate the feasibility of such optimization. The parameters were validated by a simulation using J 2 D. Based on the above results, the deformation theory is suitable for accurately predicting the occurrence of wrinkles and utilizing to optimize the influence parameters.