Finite Element Analysis of Three-Dimensional Hot Bending and Direct Quench Process Considering Phase Transformation and Temperature Distribution by Induction Heating

Abstract

In this study, a coupled thermo-mechanical-metallurgical finite element analysis (FEA) method was developed to investigate the deformation behavior in the three-dimensional hot bending and direct quench processes. In the developed FEA procedure, the temperature distribution was calculated by two methods. First was a three dimensional electromagnetic and heat conduction analysis that considered a non-linearity of permeability and magnetic transformation. Second was a simplified method that used an original heat source model for induction heating. In the deformation analysis, temperature, micro structure and strain rate dependencies of flow stress were taken into consideration. As for the microstructure evolution, an experimental formula was used to track the ferrite-austenite transformation, and Koistinen-Marburger relationship was employed to describe the austenite-martensite change. To confirm the effectiveness of the developed FEA method, the thickness change upon bending and the camber by inhomogeneous cooling were simulated. The results were in good agreement with the experimental measurements.

1. Introduction

In recent years, the automotive industry has been focusing on two issues: the development of lighter vehicles to improve fuel economy in an effort to prevent global warming, and the improvement in crash safety. A three-dimensional hot bending and direct quench (3DQ) mass processing technology of tube has been developed^1,2) as a means of satisfying these two conflicting needs. The technology enables automotive parts to be a hollow tubular structure with an ultra-high-tensile strength and a three-dimensional complex shape.

Because a high dimensional accuracy is usually required for automotive parts, it is necessary to know accurately the deformation behavior of tubes in 3DQ. However the actual 3DQ is a complicated process where electromagnetic fields, heat conduction and heat transfer, latent heat of transformation, microstructure change, mechanical properties change, thermal and transformation expansion, and transformation plasticity are involved.

In the past, some numerical models for tube bending processes using induction heating (Induction bending) have been proposed. Asao et al.³⁾ and Kuriyama et al.⁴⁾ predicted the bending moment and wall thickness of cylindrical tubes using elementary analyses. The models were based on the assumption of plane strain, and they stated that the cross-section shape did not change during the deformation.

A finite element method (FEM) has also been used to simulate the induction bending process. Wang et al.,⁵⁾ Hu et al.⁶⁾ and Tropp et al.⁷⁾ calculated thickness change and flattening by a coupled temperature-deformation FEM. Additionally, Lee et al.⁸⁾ tried to calculate a residual stress distribution. However, the above studies did not take into account phase transformation and an accurate temperature distribution generated by induction heating. Hence, these models are not enough to apply to the 3DQ process.

To generate the temperature distribution for deformation analysis, heat source model of induction heating is necessary. In previous studies, some heat source models^9,10) for induction heating were proposed. However, these models cannot be used in the present study because they are for spot heating of plate.

In this study, a reliable simulation framework for 3DQ is developed in which an accurate temperature distribution generated by induction heating and phase transformation are taken into account. In the framework, also a simplified method to calculate temperature distribution including a heat source model for local induction heating of tube is proposed to save the computation time. Finally the simulation accuracy is examined by comparing with experimental results.

2. Outline of the 3DQ Process

The characteristics of the 3DQ process and its products are explained in the literatures.^1,2) The outline is as follows. 3DQ is a consecutive forming technique that allows three-dimensional complex hollow bending and quenching at the same time. Figure 1 shows an example of a 3DQ machine. The machine consists of five components: an axial feeding device, guide rolls, an induction heater, a cooling device, and a bending device (Robot). 0.21 mass% carbon-boron steel tubes are used as the material. First, the tube is fed in the axial direction and heated to higher than Ac₃ temperature locally. Just after heating, the tube is bent by an industrial robot. At the same time, the tube is quenched to room temperature by the cooling device. To simulate the 3DQ process accurately using FEM, it is important to consider the temperature history and the metal structure change.

Fig. 1.

Schematic illustration of 3DQ machine using industrial robot. (Online version in color.)

3. Modeling of 3DQ

We propose a simulation framework for 3DQ including coupled thermo-mechanical-metallurgical FEM in this study. Figure 2 shows the procedure of the framework. This procedure consists of the following two stages.

Fig. 2.

Procedure of 3DQ analysis.

At the first stage, a temperature distribution is calculated. In this study, two methods are employed. One is the three-dimensional electromagnetic and heat conduction analysis¹⁴⁾ (Section 3.1.1). We call this model the three-dimensional model in this paper. This model is used in the case where the detail consideration of the effect of three dimensional distribution is needed (Section 4.3). Another is a simplified method that includes original heat source model and one dimensional heat conduction analysis (Section 3.1.2). This method is called the simplified method in this paper. When the temperature distribution is already known by an experiment or the three dimensional model, the simplified method is used to save time to prepare simulation models and the computation time (Sections 4.1 and 4.2).

At the second stage, the temperature distribution calculated in the first stage is assigned to integration points in finite elements. Then, the deformation analysis is carried out considering temperature, micro structure and strain rate dependencies of flow stress using the coupled thermo-mechanical-metallurgical FEM. The static implicit method in ABAQUS Standard is used in this stage.

To take into account the microstructure change during deformation, the Kunitake’s equations¹⁶⁾ as a function of material composition are used to describe the transformation temperatures, the experimental formula¹¹⁾ is used to track the ferrite-austenite transformation, and a Koistinen-Marburger relationship¹²⁾ is employed to describe the austenite-martensite change. In order to facilitate the model when the material composition is changed, not the dilatometry test result specific to the material but the above-mentioned generalized equations are used in the present study. Then, the flow stress is determined using these transformation models and a linear mixture rule.¹³⁾ These models were implemented into ABAQUS through user subroutines.

The above mentioned models are described in detail in the following sections.

3.1. Temperature Analysis

3.1.1. Three-dimensional Model

In this analysis, coupled electromagnetic and heat conduction analysis is carried out. A joule heat density due to induction heating is calculated by an electromagnetic analysis, and a temperature distribution is calculated by a heat conduction analysis. It is difficult to obtain an accurate three dimensional result if the actual quenching process is modeled as the tube is very long and moves in the longitudinal direction. To solve such a problem, the modeling area is narrowed down to locally near the induction coil. Further, the feeding of material is considered as the flow field using an advective term in the equation of heat conduction. These analysis algorithms are described in the literature.¹⁴⁾ This analysis considers the non-linearity of permeability, the magnetic transformation and the temperature dependent of thermal and electric property.

When the temperature distribution reaches a steady state, it is assigned to integration points in finite elements used in the deformation analysis.

3.1.2. Simplified Method

In the simplified method, to generate the temperature distribution, the heat conduction analysis using the heat source model is carried out.

The heat source model is based on a simple theory and an empirical equation. An electric current density of high-frequency induction heating decreases exponentially with a distance from outer surface as follows¹⁵⁾   

$$i=i_{0}exp(-x/\delta ),$$(1)

where i is an electric current density, i₀ is an electric current density of outer surface, x is a distance from outer surface, and δ is a skin depth. The distribution of the electric current density calculated from Eq. (1) is schematically shown in Fig. 3(a). When δ is small, the electric current concentrates in the outer surface of a tube. Next, a longitudinal distribution of outer surface electric current density is approximated by a normal distribution function as follows   

$$i_0=exp\left[-\frac{{\left(z-m\right)}^2}{(2/9)a^2}\right],$$(2)

where m is a parameter for adjusting the position of the peak, a is a half width of heating area in the longitudinal direction, and z is a longitudinal position. The longitudinal distribution of electric current density calculated from Eq. (2) is schematically shown in Fig. 3(b). The value of i₀ in Eq. (2) is normalized to a range between 0 and 1.

Fig. 3.

Electric current distributions in thickness and longitudinal directions.

Here, the relationship between joule heat density and electric current density is given in the form   

$$Q=i^{2}R,$$(3)

where R is an electric resistance and Q is a joule heat density. From Eqs. (1), (2), (3), a two-dimensional joule heat density is approximated by the following equation   

$$Q=Q_{max}{\left[exp\left(-\frac{{\left(z-m\right)}^2}{\left(2/9\right)a^2}\right)exp\left(-\frac{x}{\delta }\right)\right]}^2,$$(4)

where Q_max is a maximum joule heat density. The distribution of the joule heat density calculated from Eq. (4) is schematically shown in Fig. 4, where the origins of x-axis and z-axis are the outer surface of tube and the point of center of induction coil, respectively. τ is a thickness of tube.

Fig. 4.

Heat source model for local induction heating of tube. (Online version in color.)

Using the joule heat density distribution of Eq. (4), one dimensional heat conduction analysis in thickness direction is then carried out. The governing equation of one dimensional heat conduction is of the form   

$$\frac{\partial T}{\partial t}=\frac{\lambda }{\rho c}\frac{{\partial }^{2}T}{\partial x^2}+\frac{Q}{\rho c},$$(5)

where T is a temperature, t is a time, λ is a heat conduction coefficient, ρ is a density, and c is a specific heat. As for the boundary condition, convection heat transfer by water jet on the outer surface of tube is considered by   

$$-\lambda {\frac{\partial T}{\partial x}|}_{x=0}=h\left(T_s\left(t\right)-T_w\right),$$(6)

where h is a heat transfer coefficient, T_s is an outer surface temperature and T_w is a water temperature. On the inner side of tube, thermal heat insulation is assumed by the following equation   

$${\frac{\partial T}{\partial x}|}_{x=\tau }=0.$$(7)

From Eqs. (4), (5), (6), (7) the heat conduction analysis is carried out using a finite difference method with Crank-Nicolson method. Figure 5 shows the schematic illustration of finite difference analysis. The finite difference lattice is moved along the z-axis with the tube feeding speed V. Then the joule heat density Q is given by Eq. (4). As a result, a two-dimensional steady state temperature distribution in the thickness and longitudinal directions is obtained.

Fig. 5.

One dimensional heat conduction analysis using moving finite difference lattice.

The advantages of the simplified method are as follows: it is not necessary to model the coil and tube geometries, the computation time is very short and an experimental temperature distribution can be reproduced easily by adjusting the parameters in Eq. (4).

3.1.3. Comparison between Three-dimensional Model and Simplified Method

The effectiveness of simplified method is examined by the comparison with the three-dimensional model.

The analysis using the three dimensional model was carried out considering a circular cross-section tube with 1.6 mm in thickness and 52.56 mm in outer diameter. The feeding speed was 80 mm/s. The center of coil was z=0 mm and the cooling area was z>25.33 mm. The frequency of induction heating was 10 kHz. The magnetic transformation temperature (737°C) was also taken into account. The distribution of joule heat density obtained by this analysis is shown in Fig. 6. The skin depth rapidly increased at around z=0. Figure 7 shows the temperature distribution in the longitudinal direction on the inner and outer surfaces. The temperature exceeds the magnetic transformation temperature (737°C) at around z=0. Then the skin depth rapidly increased here.

Fig. 6.

Joule heat density by three dimensional model. (Online version in color.)

Fig. 7.

Comparison between three-dimensional model and simplified method.

In the simplified method, the thickness, feeding speed and cooling area were the same as those used in the three-dimensional model. The parameters of Eq. (4) were as follows: m=–5.67 mm, a=53 mm, Q_max=2.39×10¹⁰ W/m³, and δ =2.8 mm. The thermal properties were assumed constant and the values at 600°C were used: λ=35.6 W/mK, ρ= 7830 kg/m³ and c=731 J/kgK. The distribution obtained by Eq. (4) is shown in Fig. 8. The complex joule heat density distribution that was predicted in the three-dimensional model (Fig. 6) could not be achieved in the simplified method. However, the simplified method could describe the bell shape distribution in the longitudinal direction and the skin effect in the distribution in thickness direction as those predicted using the three-dimensional model.

Fig. 8.

Joule heat density by heat source model (Eq. (4)). (Online version in color.)

As shown in Fig. 7, the temperature distributions predicted using the simplified method are also consistent with those of the three-dimensional model with a sufficient accuracy for practical use. There is a little discrepancy at a temperature of around 200°C or less, but this is not a practical problem since it is a low temperature range with small deformation. The above results confirmed the effectiveness of the simplified method.

3.2. Modeling of Microstructural Evolution

For the analysis of the phase transformations, a mathematical model based on physical and empirical approaches is used. The 3DQ deformation model has both a heating area and a cooling area. In the heating area, an above-Ac₃ temperature is mandatory to obtain the full austenite structure. In the cooling area, the full martensite structure is obtained without bainite or ferrite transformation as the cooling rate in 3DQ is very high.

The transformation temperatures of Ac₁°C, Ac₃°C and M_s°C are described using Kunitake’s equations (Eqs. (8), (9), (10)¹⁶⁾) and a regression formula for chemical composition to provide a versatile model. X_i is a weight percentage of a composition i   

$$\begin{array}{l}A{c}_1=\text{727}\text{.0}-\text{32}\text{.7}\times X_{\text{C}}+\text{14}\text{.9}\times X_{\text{Si}}+\text{2}\text{.0}\times X_{\text{Mn}}\\ -\text{17}\text{.0}\times X_{\text{Cu}}-\text{14}\text{.2}\times X_{\text{Ni}}+\text{17}\text{.8}\times X_{\text{Cr}}+\text{25}\text{.6}\times X_{\text{Mo}},\end{array}$$(8)

  

$$\begin{array}{l}A{c}_3=\text{912}\text{.0}-\text{230}\text{.5}\times X_{\text{C}}+\text{31}\text{.6}\times X_{\text{Si}}-\text{20}\text{.4}\times X_{\text{Mn}}\\ -\text{39}\text{.8}\times X_{\text{Cu}}-\text{18}\text{.1}\times X_{\text{Ni}}-\text{14}\text{.8}\times X_{\text{Cr}}+\text{16}\text{.8}\times X_{\text{Mo}},\end{array}$$(9)

  

$$\begin{array}{l}{M}_S=\text{560}\text{.5}-\text{407}\text{.3}\times X_{\text{C}}-\text{7}\text{.3}\times X_{\text{Si}}-\text{37}\text{.8}\times X_{\text{Mn}}\\ -\text{20}\text{.5}\times X_{\text{Cu}}-\text{19}\text{.5}\times X_{\text{Ni}}-\text{19}\text{.8}\times X_{\text{Cr}}-\text{4}\text{.5}\times X_{\text{Mo}}.\end{array}$$(10)

The approximation Eq. (11) is used to simulate the formation of austenite in the heating area¹¹⁾   

$${\xi }_{\gamma }=1-exp\left[-cA{\left(\frac{T-A{c}_1}{A{c}_3-A{c}_1}\right)}^{da}\right]\hspace{1em}(A{c}_1<T<A{c}_3),$$(11)

where ξ_γ is a volume fraction of austenite and cA and da are constants.

At the cooling area, the transformation from austenite to martensite is described by Koistinen-Marburger’s formula¹²⁾   

$${\xi }_M=1-exp\left[-0.011\left(M_S-T\right)\right]\hspace{1em}(T\le M_S),$$(12)

where ξ_M is a volume fraction of martensite. An austenite volume fraction ξ_γ in the cooling area is defined as follows.   

$${\xi }_{\gamma }=1-{\xi }_M\hspace{1em}(T\le M_s).$$(13)

3.3. Temperature, Metal-structure, and Strain Rate Dependencies of Flow Stress

To carry out an elastic-plastic deformation analysis considering temperature, metal-structure and strain rate dependencies of flow stress, the flow stresses were measured experimentally at various conditions under uniaxial tension using electrical resistance heating and inert gas flow cooling as follows.

To measure the flow stress of ferrite-pearlite structure, the temperature was increased by heating rate of 500°C/s. When the temperature reached the testing temperature, the tensile test was immediately carried out. The tests were carried out in the temperature range between room temperature and 800°C.

To measure the flow stress of austenite structure, the temperature was first increased to over 1000°C to obtain full austenite structure. After that, the specimen was cooled down in 100°C/s until the temperature reached a preset value and the tensile test was immediately carried out. The upper critical cooling rate of this material is 30°C/s;¹⁷⁾ hence the tensile test can be carried out in full supercooled austenite state. The test was conducted in the temperature range between 400 and 1100°C.

To measure the flow stress of the martensite structure, quenched full martensite material was used. The tensile test was carried out after the temperature was increased to a preset value. The test was conducted in the temperature range between room temperature and 400°C.

In this deformation analysis, mechanical properties such as Young’s modulus and flow stress of mixture are defined by the linear mixture rule.¹³⁾ For example, the flow stress of mixture is given as follows   

$${\sigma }_y\left(T\right)=\sum _{I=1}^N{\sigma }_{yI}\left(T\right){\xi }_I,$$(14)

where σ_yΙ(T) is an experimental value of flow stress of phase I at the temperature, ξ_Ι is a predicted value of the volume fraction of the phase, and σ_y(T) is a calculated flow stress. In this analysis, N is 3 because the three kinds of phases appear.

Strain rate dependency of flow stress is defined by the following method. The tensile test was carried out in the conditions: temperature of 800°C to 1000°C and strain rate of less than 0.6 s^–1.

Figure 9 shows the obtained flow stresses at a strain of 0.03: the medium value of this simulation. Here, the experimental result is approximated by the following Cowper-Symonds law¹⁸⁾   

$$\sigma ={\sigma }_0\left\{1+{\left(\frac{\dot{\epsilon }}{D}\right)}^{1/P}\right\},$$(15)

where D and P are material constants, $\dot{\epsilon }$ is a strain rate, and σ₀ is a flow stress at $\dot{\epsilon }$ = 0. The lines in Fig. 9 are the approximations using Eq. (15). In this case, the increase of flow stress from 800°C to 1000°C is described by a pair of parameters: D=1.53 s^–1, and P=2.57.

Fig. 9.

Flow stress fitting by Cowper-Symonds law.

3.4. Thermal and Transformation Strain

Thermal and transformation strain is defined as follows. First, a density of each metal structure is calculated by Miettinen’s equations (Eqs. (16), (17), (18), (19), (20), (21), (22)¹⁹⁾).   

$$\begin{array}{l}{\rho }^{\gamma }\left(T\right)={\rho }_{Fe}^{\gamma }+\left(-118.26+7.39\times {10}^{-3}T\right)X_C-68.24X_{Si}\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }-6.01X_{Mn}+12.45X_{Mo}+(-7.59+3.422\times {10}^{-3}T\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }-5.388\times {10}^{-7}{T}^2-1.4271\times {10}^{-2}{X}_{Cr})X_{Cr}\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }+(1.54+2.267\times {10}^{-3}T-11.26\times {10}^{-7}{T}^2\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }+6.2642\times {10}^{-2}{X}_{Ni})X_{Ni},\end{array}$$(16)

  

$${\rho }^{\alpha +C}(T)={\left(\frac{6.69-X_C}{6.69{\rho }_{Fe}^{\alpha }\left(T\right)}+\frac{X_C}{6.69{\rho }^C\left(T\right)}\right)}^{-1}+{\rho }_{sub}\left(T\right),$$(17)

  

$${\rho }^M(T)={\rho }_{Fe}^{\alpha }\left(T\right)+\left(7.92X_C-168.2\right)X_C+{\rho }_{sub}\left(T\right),$$(18)

  

$${\rho }_{Fe}^{\gamma }(T)=8\mathrm{\hspace{0.17em} }099.79-0.506T,$$(19)

  

$${\rho }_{Fe}^{\alpha }(T)=7\mathrm{\hspace{0.17em} }875.96-0.297T-5.62\times {10}^{-5}{T}^2,$$(20)

  

$${\rho }^C(T)=7\mathrm{\hspace{0.17em} }686.45-0.0663T-3.12\times {10}^{-4}{T}^2,$$(21)

  

$$\begin{array}{l}{\rho }_{sub}\left(T\right)=-36.86X_{Si}-7.24X_{Mn}+30.78X_{Mo}\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }+(-8.58+1.229\times {10}^{-3}T+8.52\times {10}^{-8}{T}^2\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }+1.8367\times {10}^{-2}{X}_{Cr})X_{Cr}+(-0.22-4.7\times {10}^{-4}T\\ \hspace{1em}\hspace{1em}\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }\mathrm{\hspace{0.17em} }-1.855\times {10}^{-7}{T}^2+0.104608X_{Ni})X_{Ni},\end{array}$$(22)

where ρ^γ(T) is a density of austenite structure, ρ^α+C(T) is a density of ferrite-pearlite structure and ρ^M(T) is a density of martensite structure. ${\rho }_{Fe}^{\gamma }$ (T) is a density of austenite structure of pure iron, ${\rho }_{Fe}^{\alpha }$ (T) is a density of ferrite structure of pure iron and ρ^C(T) is a density of cementite structure. ρ_sub(T) is a member to express the effect of composition except carbon.

Using the above definitions, the thermal and transformation strain dε_T,Tr is given in the form¹¹⁾   

$$d{\epsilon }_{T,Tr}=\frac{l^{k+1}\left(T\right)}{l^k\left(T\right)}-1=\sqrt[3]{\frac{{\rho }^k\left(T\right)}{{\rho }^{k+1}\left(T\right)}}-1,$$(23)

where l(T) is a length of line element in the temperature T. ρ ^k(T) and ρ ^κ+1(T) are densities at the beginning and end of increment, respectively. ρ ^κ+1(T) is calculated using the linear mixture rule with volume fraction ξ_i (T) and the density ρ_i (T) at the end of increment as follows.   

$${\rho }^{k+1}\left(T\right)=\sum _{i=1}^n\left\{{\xi }_i\left(T\right){\rho }_i\left(T\right)\right\}.$$(24)

Figure 10 shows the temperature-strain curves obtained from the simulation using Eqs. (8), (9), (10), (11), (12), (13) and Eqs. (16), (17), (18), (19), (20), (21), (22), (23), (24) and the experiment. In the experimental dilatometry test, the heating rate of 10°C/s and the cooling rate of 60°C/s were employed. An external force was not applied. Although the simulated result is in agreement with the experimental one qualitatively, this does not match quantitatively. Especially, a discrepancy is observed in the transformation curve from austenite to martensite. A discrepancy is also observed in the slope of the cooling curve. The effect of the discrepancies on the simulation of 3DQ will be discussed in sections 4.1 and 4.2.

Fig. 10.

Temperature-strain curves obtained from simulation and experiment.

4. Results

4.1. Thickness

To confirm the analysis accuracy, an experiment and an analysis of bending processes were carried out under the three conditions shown in Table 1. Because of the validation of large deformation, the thermal and transformation strains are neglected in these analyses. The temperature distribution was calculated by the simplified method to save time for the preparation of simulation models. The parameters in heat source model in section 3.1.2 were identified as m=–5.67 mm, a=53 mm and δ=2.8 mm. Q_max was adjusted to reach the outer surface temperature 1000°C. The following parameters are the same among the three conditions. The feeding speed is 80 mm/s, and a heat transfer coefficient of the cooling area is 25000 W/m²K. Figure 11 shows the geometry and boundary conditions of FE model. A four-node shell element with reduced integration was used. Five integration points were given through thickness. The number of element is about 20000. The movement of robot was modeled and reproduced as displacement boundary conditions given at the top of tube.

Table 1. Analytical and experimental conditions.

	Tube geometry	Bending radius
Condition A	40 H × 46 W × 1.6 t	1000 mm
Condition B	24 H × 56 W × 2.3 t	1120 mm
Condition C	24 H × 56 W × 2.3 t	250 mm

Fig. 11.

Geometry and boundary condition of FE model (Condition C).

Figure 12 shows the thickness strain distribution in the longitudinal direction at the center of inner and outer side of tubes (points a and b). The simulation results are in good agreement with the experimental ones.

Fig. 12.

Comparison of longitudinal thickness strain distributions between experiment and FEA (Condition C).

Figure 13 shows the comparison of maximum and minimum thickness strains between the experiment and the simulation for all the conditions. Clearly the simulation results are in good agreement with the experiment results regardless of the bending conditions.

Fig. 13.

Comparison of thickness strain between experiment and FEA.

In these analyses, the plastic deformation occurs in the portions of the ferrite-pearlite structure, the austenite structure and the mixture of ferrite-pearlite and austenite structure. The transformation model describes well the transformation from ferrite-pearlite structure to austenite structure as shown in Fig. 10, leading to the accurate prediction of flow stress change. On the other hand, because the plastic deformation does not occur after the martensite transformation starts, it is also considered that the martensite transformation in which the discrepancy was observed in Fig. 10 does not affect the simulation results.

4.2. Camber by Inhomogeneous Cooling

The above bending analysis was good enough to confirm the analysis accuracy of large deformation. However, this model was not suitable to confirm the thermal and transformation deformation because the deformation due to bending is much larger than thermal deformation. Therefore, quenching analyses with various inhomogeneous cooling conditions were carried out to confirm the accuracy of the simplified method and the thermal and transformation strain models.

Figures 14 and 15 show the boundary conditions and geometry of the FE model, respectively. They were determined based on practical 3DQ conditions.

Fig. 14.

Geometry and boundary condition of FE model (Temperature distribution for cooling condition 3).

Fig. 15.

Top view of tube.

A rectangular tube of 46 mm in height, 40 mm in width, 1.6 mm in thickness, and 1000 mm in length was used. Feeding speed was 80 mm/s.

Table 2 shows the simulation conditions. To clarify the transformation effect, the analysis was carried out with or without considering the transformation. When the transformation was not considered, metal structure remained in the ferrite-pearlite structure (Cases 7 to 12 in Table 2). In these cases, the flow stress of ferrite-pearlite structure at 800°C or over was substituted by the flow stress of austenite because the flow stress of ferrite-pearlite structure is not be able to obtain at 800°C or over.

Table 2. Simulation conditions.

Case No.	Transformation	Constraint condition	Cooling condition
1	Yes	Robot	1
2	Yes	Robot	2
3	Yes	Robot	3
4	Yes	Free	1
5	Yes	Free	2
6	Yes	Free	3
7	No	Robot	1
8	No	Robot	2
9	No	Robot	3
10	No	Free	1
11	No	Free	2
12	No	Free	3

To examine the analysis accuracy under different boundary conditions two kinds of displacement boundary conditions were employed as shown in Fig. 14 and Table 2, i.e. the edge of tube was constrained by the robot or was free.

When the front edge was constrained by the robot motion (Cases 1 to 3 and 7 to 9 in Table 2), the longitudinal displacement Ux was set to free and the other degrees of freedom were fixed. On the other hand, in Cases 4 to 6 and 10 to 12, all degrees of freedom were set to be free.

To change the curvature of camber, three kinds of cooling conditions were employed.

Condition 1 is a basic condition where homogeneous cooling is assumed. In conditions 2 and 3, inhomogeneous cooling is assumed in which the cooling rate in plane A (Fig. 15) is lower than the other planes. The cooling rate in plane A of condition 3 is lower than condition 2.

In the experiment, the amount of water flow onto plane A was reduced to decrease the cooling rate in condition 2 and 3. It should be noted that this inhomogeneous cooling condition is only for the verification analysis and is not in practical use of 3DQ. In the analyses, a cooling curve was fitted to experimental data by try-and-error adjustment using the simplified method. In such cases, the simplified method is effective to save the computation time.

Figure 16 shows the temperature evolutions at the inner surface of the tube in planes A, B, C and D (Fig. 15) obtained by the experiment and the simulation. The data was evaluated along the lines L_A and L_C (Fig. 14). The simulation results agree well with the experimental results.

Fig. 16.

Temperature evolution of each cooling condition.

Curvature after cooling was evaluated by three-point curvature in a length of 800 mm.

Figure 17(a) shows the experimental results of curvature after quenching. A positive curvature means plane A is the inner side of camber and a negative curvature the outer side. These results have the following features. The absolute curvature becomes larger when the cooling rate of plane A becomes lower. And the curvature variation with the robot constraint condition is opposite from those of the free condition. The absolute values of curvature obtained with the robot constraint conditions are smaller than those of free condition.

Fig. 17.

Experimental and simulated results of camber for inhomogeneous cooling.

Figure 17(b) shows the simulated results in which phase transformation was taken into account (Cases 1 to 6 in Table 2). These results are in good agreement with the experiment.

On the other hand, the simulated results in which phase transformation was not taken into account (Fig. 17(c)) do not agree with the experimental results. In the free condition (Cases 10 to 12 Table 2), the direction of camber was opposite from that observed in Fig. 17(b). In the robot constraint condition (Case 7 to 9 of Table 2), the curvature became much larger than that observed in Fig. 17(b).

The above results confirm the effectiveness of this analysis framework. Furthermore, it turned out that an accurate thermal deformation analysis is necessary by considering a temperature distribution, a metal-structure change and a mechanical boundary condition.

Because the results shown in Fig. 17(b) are in good agreement with the experiment, the effect of the discrepancies in the temperature-strain curve (Fig. 10) on the simulation results would be small as far as we tested in the present study. However, there is a room to improve the accuracy of the transformation model. For example, Bok et al.²⁰⁾ proposed a model which can approximate the volume fraction of martensite using some fitting parameters. Using this model, it is expected that the simulation accuracy with the variety of forming conditions may be improved. This will be our future work.

4.3. Deformation Behavior in 3DQ

To clarify the deformation and the stress evolution during 3DQ, a simulation of a simple bending process was carried out. First, the features of temperature distribution generated for the bending analysis are described in section 4.3.1. After that the bending analysis conditions and the results are described in section 4.3.2.

4.3.1. Temperature Distribution for Bending Analysis

To generate the three dimensional temperature distribution for the bending analysis, the three-dimensional model described in section 3.1.1 was used. The analysis was carried out under the following conditions: a rectangular tube of 35 mm in height, 45 mm in width, 1.6 mm in thickness, 80 mm/s in feeding speed, and maximum temperature of 1000°C. The heat transfer coefficient of the cooling area was assumed to be 20000 W/m²K. Figures 18 and 19 show the results of joule heat density and temperature distributions, respectively. The broken lines show the coil portion and cooling area. The peak of joule heat density appears at the upper stream of center of coil (point A in Fig. 18). This is because of the decreasing of heat efficiency owing to the magnetic transformation which occurs at 737°C during heating. On the other hand, the temperature is distributed homogeneously in the circumferential direction. The peak temperature arises just before the cooling area (point B in Fig. 19).

Fig. 18.

Joule heat density at outer surface obtained by three dimensional model.

Fig. 19.

Temperature distribution at outer surface obtained by three dimensional model.

4.3.2. Deformation Analysis

Figure 20 shows the geometry and boundary conditions of FE model. The length of the tube was 1000 mm, and the bending radius was 1000 mm at the center line of tube. The element size was 1.8 mm each side. To simplify the simulation, thermal and transformation strains were neglected. The temperature distribution shown in Fig. 19 was applied.

Fig. 20.

Geometry and boundary condition of FE model.

Figure 21 shows the result at the center of thickness of the evaluation nodes in Fig. 20. The longitudinal distributions of the temperature, equivalent stress, subsequent yield stress, and longitudinal strain rate are shown. The subsequent yield stress decreases rapidly from 400 N/mm² to 40 N/mm² as the temperature increases from 600°C to 800°C. The plastic deformation occurs at the portion where the temperature is higher than 800°C. There after the plastic deformation does not occur because the material is cooled. As a result, the plastic deformation occurs only within the narrow area with a longitudinal length of about 20 mm. Such a characteristic of the plastic deformation is effective to suppress the change of cross sectional shape. The yield stress is increased by martensite transformation.

Fig. 21.

Distributions of temperature, strain rate, and stress at center of thickness.

The equivalent stress is relatively small in the whole region regardless of the temperature. This shows that the magnitude of stress in the whole region is determined by the flow stress in the hot region, yielding small springback in the induction bending.

5. Conclusions

In this study, a simulation framework of the three-dimensional hot bending and direct quench (3DQ) process is proposed to investigate the deformation behavior. The features of this framework and the simulation results are as follows.

(1) A coupled thermo-mechanical-metallurgical FEM analysis framework was developed. In the electromagnetic and temperature analysis, a non-linearity of permeability, magnetic transformation and temperature dependent of thermal and electric property were considered. The three dimensional temperature distribution can be obtained by this analysis. In the deformation analysis, temperature, micro structure and strain rate dependencies of flow stress were taken into consideration. The phase transformation rate was calculated using Kunitake’s equation and Koistinen-Marburger relationship. The flow stress of the mixture was calculated by a linear mixture rule. A strain rate dependent of flow stress was defined by Cowper-Symonds law. To calculate the thermal and transformation strain, Miettinen’s density law was employed.

(2) A heat source model for local induction heating of tube and the simplified method to generate the temperature distribution was also proposed. The proposed method is effective in terms of the computation time.

(3) The analysis accuracy of large deformation was confirmed through a bending analysis. And the thermal deformation accuracy was confirmed through an analysis of quenching process with inhomogeneous cooling.

(4) The detail of deformation mechanism in 3DQ was revealed by this analysis. For example it was found that the plastic deformation area was narrow of about 20 mm, and the deformation stress was about 40 N/mm².

(5) Although the simulated results of 3DQ are in agreement with the experimental ones, the discrepancies were observed in the dilatometry test. The improvement of the thermal and transformation models will be our future work.

References

• 1)   A.  Tomizawa,  N.  Shimada,  S.  Inoue,  H.  Kikuchi and  S.  Kuwayama: Proc. 2010 Jpn. Spring Conf. Technol. Plast., JSTP, Tokyo, (2010), 207.
• 2)   A.  Tomizawa,  N.  Shimada and  H.  Matsuda,  H.  Mori,  M.  Hara and  S.  Kuwayama: Proc. Int. Conf. Hydroforming of Sheets, Tubes and Profiles, MAT INFO Werkstoff-Informationsgesellschaft mbH, Frankfurt, (2012), 127.
• 3)   H.  Asao,  K.  Okada,  H.  Yonemura and  S.  Fujishima: J. Jpn. Soc. Technol. Plast., 33 (1992), No. 372, 49.
• 4)   S.  Kuriyama and  T.  Aida: Proc. of Advanced Technology of Plasticity, Vol. 1, ICTP, Beijing, (1993), 464.
• 5)   Z.  Wang and  Z.  Hu: J. Mater. Proc. Technol., 21 (1990), 275.
• 6)   Z.  Hu and  J. Q.  Li: J. Mater. Proc. Technol., 91 (1999), 75.
• 7)   M.  Tropp,  M.  Handrlk,  P.  Kopas and  M.  Saga: Pipeline Gas J., 239 (2012), No. 11, 86.
• 8)   Y.  Lee,  J.  Lim and  Y.  Moon: Mater. Trans., 53 (2012), No. 5, 847.
• 9)   K.  Tanaka,  S.  Imatani and  T.  Inoue: Proc. Jpn. Cong. on Mater. Research, Vol. 42, Soc. Mat. Sci. Jpn., Tokyo, (1998), 251.
• 10)   Y.  Luo,  M.  Ishiyama and  H.  Murakawa: Trans. JWRI, 33 (2004), No. 1, 59.
• 11)   K.  Okamura: J. Soc. Mater. Sci. Jpn., 55 (2006), No. 5, 529.
• 12)   D. P.  Koistinen and  R. E.  Marburger: Acta Metall., 7 (1959), No. 1, 59.
• 13)   T.  Inoue,  K.  Tanaka and  S.  Nagaki: Analysis of Solid Mechanics and Phase Transformation, Taiga publishing, Tokyo, (1995), 82.
• 14)   N.  Okada,  H.  Sakamoto,  T.  Matsuo and  K.  Takatani: Tetsu-to-Hagané, 98 (2012), No. 7, 38.
• 15)   A.  Chikami: Denki-keisan, 48 (1980), No. 14, 196.
• 16)   T.  Kunitake: J. Jpn. Soc. Heat Treat., 41 (2001), No. 3, 164.
• 17)   T.  Nishibata and  N.  Kojima: Tetsu-to-Hagané, 96 (2010), No. 6, 26.
• 18)   G. R.  Cowper and  P. S.  Symonds: Brown Univ. Div. of Applied Mech. Rep., 28 (1952).
• 19)   J.  Miettinen: Metall. Mater. Trans. B, 28B (1997), 281.
• 20)   H.-H.  Bok,  J.-W.  Choi,  F.  Barlat,  D. W.  Suh and  M.-G.  Lee: Int. J. Plast., 58 (2014), 154.