Abstract
This paper describes the procedure of electromagnetic and thermal modeling of an induction hardening load of 42CrMo4 steel. In addition to the temperature and field level dependent physical properties of the material, the magnetic behavior of the load is captured by means of some kind of non-linear impedance boundary condition which simplifies the computational cost of the simulation. The numerical results show the critical behavior around Curie temperature. Finally, a comparison between several simulation results and experimental measurements are provided to assess the usefulness of the proposed electro-thermal simulation.
Introduction
Since the end of the nineteenth century, Induction Heating (IH) technology has been developed, and thus, increasingly employed [1]. It is mostly used in industrial, domestic, and medical applications, but its potential is far from limited [2].
Especially focusing on industrial applications [3], the same workpiece can be treated under different operating conditions, e.g., a higher frequency is needed for induction hardening than for melting. The operation point is determined by the current flowing through the inductor. In addition, not only the heat treatment type but the characteristics of the load also vary depending on the working frequency and current level.
Temperature is also a factor that modifies the properties of materials, particularly at temperatures above Curie point (700–800 °C for typical steels). Such change of properties directly influences the electrical behavior of the inductor-load pair, which must be predicted to properly analyze the power converter supplying the heating system. The basis of an industrial hardening process is to exceed the austenizing temperature in the area to be processed, to afterwards create the martensite layer, the hardened one, through the rapid cooling of the austenite. It is worth mentioning that this austenizing point is above the Curie temperature; for this reason, modeling the load at high temperatures is a key factor when designing and foreseeing this type of IH processing.
Electro-thermal modeling of the system permits accurate analysis of such behavior, but the computational cost, in that case, is very high and requires precise characterization of the properties of the load. Numerical computations based on Finite Element Method (FEM) are regularly employed and commonly used to model IH systems [4–6].
In this work, the 42CrMo4 steel part is modeled including the dependence of its properties on magnetic field and temperature. The equivalent magnetic permeability depends on both parameters [7], while the other electrical and thermal properties depend only on temperature [7,8]. The preceding dependencies have already been studied in several papers in the literature. This article combines the mentioned dependencies in a single FEM analysis using a new approach based on a non-linear boundary condition description. This approach allows a realistic description of the involved phenomena with a reduced computational cost.
After this introduction, the electro-thermal model under study will be explained; starting with the magnetic, electrical and thermal properties of the 42CrMo4 steel, and afterwards, detailing the FEM modeling process together with some results. In the third section, the model is validated comparing some simulation data with some experimental measurements. To finish, some conclusions will summarize the paper content.
Electro-thermal model
In this work, a typical industrial IH load geometry composed of a cylindrical workpiece heated up by a 6-turn coil, as can be seen in Fig. 1, is analyzed. The material of the load is 42CrMo4 steel, so-called AISI 4140, which is commonly used in industry to manufacture components of different machines for a wide range of applications. This steel will be the reference material in the analysis performed in this paper, including the properties provided in previous literature studies [7–9].
Induction hardening experimental setup.
In the metalworking industry, 42CrMo4 carbon steel is extensively used, but manufacturers do not usually provide complete information about its physical and mechanical properties which are often given at ambient temperature. Even though, in [7,8], the temperature-dependent parameters can be found under certain conditions of interest. An adequate characterization of the characteristics of the load is needed for a proper modeling of the hardening process due to the high temperature reached. The high-temperature level involved in this kind of processes entails relevant changes in the material properties.
In a first step, the characterization of the magnetic properties of the material at elevated temperature levels is a complex task but plays an essential role in the complete modeling of the load. Researchers from Lund University have been working in this area, gaining expertise in this kind of measurements [7,10]. The characteristics measured by this research group correspond to the temperature dependent BH curve, neglecting the hysteresis effects in order to simplify the analysis. The reference analytical dependence of the magnetic characteristics was previously proposed in [11]:
The BH curve strongly influences the load modeling due to the presence of non-linearity, in particular, associated with magnetic saturation, which implies a dependence on the coil current level.
The electromagnetic properties of the material also encompass the electrical conductivity, 𝜎, which will be expressed by its inverse, i.e., the electrical resistivity, 𝜌. As for the magnetic properties, a dependence of electrical resistivity on temperature will be considered, 𝜌(T).
Figure 2 shows the physical properties of the 42CrMo4 steel from the experimental results provided in [7]. BH curves at different temperatures are represented in Fig. 2a. The electrical resistivity of the material is linearly dependent on two different ranges, as can be seen in Fig. 2b.
Magnetic and electrical properties of the 42CrMo4 steel: (a) BH curves at different temperatures, (b) temperature dependent electrical resistivity.
Results depicted in Fig. 2b can be described by the expression given in Eq. (1). Coefficients for fitting the BH curves at different temperatures of 42CrMo4 steel are magnetic saturation, B s = 1,32 T, initial permeability, 𝜇r, i = 1860, Curie temperature, T c = 783 °C, and factor controlling the temperature dependence, C = 33,6 °C [9].
On the other hand, the electrical resistivity represented in Fig. 2b can be expressed by two linear expressions in the temperature ranges below and above 800 °C [7], which is close to the Curie temperature, as it is given as follows:
The complete electro-thermal modeling also needs to include the thermal properties of the material, namely the thermal conductivity, k, and the heat capacity, C P , whose temperature dependence is plotted in Fig. 3a and Fig. 3b, respectively [8]. In both cases, the parameters increase around Curie temperature due to the influence of the phase transition in the material, also associated with the structural changes in the hardening processing. The effect of the phase transition particularly affects the heat capacity; therefore, more power has to be provided to increase the temperature around this point.
Thermal properties of the 42CrMo4 steel: (a) thermal conductivity, (b) heat capacity.
An induction heating process is a complex mixture of three phenomena: electromagnetism, heat transfer, and metallurgy techniques, thus, multiphysics analyses are needed to achieve accurate estimations. Hence, numerical methods based on finite element analysis are commonly used to model IH systems. This work will be focused on the electro-thermal study, considering the mutual coupling in a single simulation. A commercial finite element method tool, COMSOL Multiphysics [12], is used to simulate the numerical model of the IH system.
Induction heating characteristics
Induction heating characteristics
Two-dimensional axisymmetric meshed geometry for FEM simulation: (a) electromagnetic domain, (b) thermal domain.
The summary of the IH system characteristics is listed in Table 1. The load is a cylindrical billet of 42CrMo4 carbon steel, and the coil is a copper hollow tube of 6 turns. The model is simulated in a two-dimensional axisymmetric environment as shown in Fig. 4.
Convection and radiation power losses of the load are modeled including the billet convection coefficient, h = 10 W/m2K [13,14], and emissivity, 𝜖 = 0,4 [15]. However, due to the high-power density levels transferred to the load by the induction process, the power losses of the two preceding phenomena have a moderate impact on the final results.
The electromagnetic modeling is performed in the frequency domain, and the thermal modeling is done in the time domain to capture the transient behavior of the load. Both types of simulations can be coupled using two COMSOL modules. However, the electromagnetic modeling of the load is a very complex task due to the BH curve saturation, but to facilitate the processing, a non-linear impedance boundary condition has been adopted. This condition allows to replace a domain by a boundary with special properties. Therefore, the mesh inside the billet is not required for the electromagnetic analysis and the computational cost is simplified (Fig. 4a). However, a mapped mesh is used for the thermal analysis (Fig. 4b). It is worth mentioning that the thermal domain mesh could be coarser, but a finer mesh have been employed to gain precision, as the computational cost is not highly elevated.
The equivalent magnetic permeability of the load at a certain temperature has been obtained by numerically solving the diffusion equation of the magnetic field for a BH curve given by Eq. (1), a similar methodology of the one proposed in [16,17]. The temperature dependence is included in the equivalent impedance by a factor similar to the exponential shown in Eq. (1). It should be noted that the preceding procedure provides a magnetic field and temperature dependent complex relative magnetic permeability, 𝜇 r (H, T) which can be inserted in the COMSOL simulation, as a physical property of the load.
The simulation is performed for a transient of 20 s, feeding the coil by a single-harmonic current of 12,5 kHz and around 600 A of amplitude. The electromagnetic properties are updated when the temperature distribution changes, then, some intermediate frequency-domain simulations are done in selected instants by the FEA tool.
As a result, the equivalent impedance of the inductor-load system depends on the time. Figure 5a and Fig. 5b show the transient dependence of the equivalent resistance and inductance, respectively, at the working frequency. It can be observed that around t = 7 s the equivalent impedance exhibits a rapid change because the Curie temperature is reached in the area located in the middle of the billet, therefore, a local saturation of the material can imply a considerable effect in the electrical characteristics of the IH load. The preceding behavior needs to be considered in the power electronics control in order to ensure that the dynamics of the close-loop current control is capable to handle the change on the electrical properties of the load.
Equivalent impedance of the IH system: (a) resistance, R eq , (b) inductance, L eq .
Figure 6a shows that even stabilizing the value of the current supplied to the inductor, the change in the equivalent resistance implies a time dependence on the power transferred to the load. At the same time, some power losses associated with convection and radiation slightly reduce the net energy accumulated by the load, as can be seen in Fig. 6b, but the temperature distribution in the load is mainly determined by the surface power density transferred to the load. The simulation results for the power losses in the load are based on constant convective and emissivity coefficients, but accurate results can be achieved by including temperature dependence.
Power flux in the load: (a) power induced into the load, (b) convection and radiation power losses in the billet.
Temperature distribution (°C) in the billet at t = 1, 5, 10, and 15 s.
The temperature distribution in the load is represented in Fig. 7 for some selected times. Higher temperatures are initially reached close to the surface of the middle area of the billet. Next, the deep regions are also around Curie temperature. Temperatures significantly above the Curie point are not observed due to magnetic saturation as well as the increase in heat capacity.
The heat flux in the load is represented in Fig. 8. As can be observed, at surface temperatures below Curie temperature, the power flux reaches the maximum values, but after the magnetic saturation of the material, the performance of the power transfer by induction clearly reduces. Therefore, the phase transition associated with the magnetic saturation contributes to the homogenization of the temperature profile of the billet.
In a first step, a small-signal measurement of the impedance of the inductor with a 42CrMo4 carbon steel load has been carried out. The measurement has been performed with a high-precision Agilent 4284A LCR meter at frequencies from 20 Hz to 100 kHz. A comparison with the simulation results from the COMSOL model is depicted in Fig. 9a and Fig. 9b for the resistance and inductance, respectively. As can be seen, the agreement between experimental measurements and simulation data is quite good, validating the geometrical modeling of the induction heating system.
Heat flux (W/mm2) in the workpiece at t = 1, 5, 10, and 15 s.
Comparison of simulation and small signal measured equivalent impedance of the IH system: (a) resistance, R eq , (b) inductance, L eq .
Next, an experiment at rated power level is conducted. The temperature of the middle point of the billet is measured by means of a K-type thermocouple for a transient of 20 seconds. The experimental temperature results are compared to simulation results in Fig. 10.
Comparison of simulation and experimental measured temperature at the boundary middle point of the billet of the IH system.
The surface central point of the billet in the simulation reaches the Curie point at around t = 7 s, and then the temperature almost becomes stable. The reason to explain this behavior is twofold: the material at this point saturates, as it is shown in the BH curve, thus, the surface power density decreases in the close area to the middle point, and, at this temperature, the phase transition implies an increase of the heat capacity which requires a higher power flux to increase the temperature.
Although the measured temperature appears to differ from the simulated results, a couple of facts can explain this error. On the one hand, the thermocouple could affect the electromagnetic field in the close region and therefore the power density transferred to the load, the saturation at the measured point, and also induce some interference in the sensor. On the other hand, a small offset in the position modifies the time dependence of the temperature at the measured point. Furthermore, considering that thermocouples exhibit a degree of inertia, they may introduce a time lag in rapidly changing temperature situations, which frequently happens during IH processes.
In this paper, an electro-thermal simulation model of processing of a 42CrMo4 steel billet including dependence of the physical properties on the temperature as well as field levels has been analyzed and experimentally tested. An accurate characterization of the properties is critical to achieving good results. Moreover, the approach of substituting the electrical and magnetic properties of the load by an equivalent impedance boundary condition dependent on the temperature as well as the field level has proven to be a useful tool to achieve reliable simulation results.