# Mathematical Modeling of Eddy-Current Loss for a New Induction Heating Device

Authors: Hai Du; Junyuan Li; Yanbin Qu
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

A new induction heating device is presented in this paper. This device can convert mechanical energy into heat energy by utilizing eddy currents, which are induced by rotating permanent magnets. A mathematical model is established for estimating eddy-current loss of the device. The distribution of induced currents and the resultant magnetic field intensity are considered in the process of modeling the eddy-current loss and so is the mutual influence of the electric field between neighborhood pole projection areas. Particularly, the skin effect is considered by correcting the numerical integral domain of eddy current density, which has great effect on the calculating results. Based on specific examples, the effectiveness and correctness of proposed model are proved by finite element analysis. The results show that the mathematical model can provide important reference for design and structure optimization of the device.


---

## Body

## 1. Introduction

In the process of induced heating, the thermal power is due to eddy currents, when conductor is subjected to an alternating-current (AC) magnetic field. Most induction heating devices, with electric energy as input, produce the AC magnetic field by using coils [1–5]. In this paper, we present a novel induction heating device, with rotating mechanical energy as input, and the eddy currents are induced by the relative motion between conductor and permanent magnets.This device is mainly composed of a stator and a dual-rotor. The stator is made of copper, inside which there are heat-exchange pipes with an inlet and an outlet, and the dual-rotor is made of iron. The structure of the device is shown in Figure1. There is an alternative magnetic field between two rotors produced by permanent magnets with alternating pole configurations. When the dual-rotor and hence the magnetic field rotate, the fluid flowing through the heat-exchange pipes can be heated by the eddy-current loss of the stator. During this process, input rotating mechanical energy is converted into thermal energy. This device is similar to the eddy-current brake in structure and energy conversion. Many efforts have been done to the modeling of the eddy-current brake so far [6–8].Structure of the heating device. (a) Stator and dual-rotor and (b) heat-exchange pipes of stator.
(a)
(b)Based on the eddy-current brake concept, we present a mathematical model of eddy-current loss for the induction heating device at low-speed operation. First, the cylindrical induction heating device with large curvature radius is simplified as linear model, so rotating movement is changed into linear movement. Second, the surface charge, due to the motional induced electric field, at the edge of the pole projection area is obtained. From these charges, the Coulomb field intensity is calculated by using Coulomb’s law. Meanwhile, the mutual influence between neighborhood Coulomb fields is considered. Third, the eddy current in the pole projection area, due to rotational magnetic field of the permanent magnets, is derived. Finally, the brake force is calculated by applying the Lorentz force law through numerical integral method. Based on the relationship between the brake force and the loss power, the eddy-current loss is obtained at different speeds. Since the eddy current generates the magnetic flux opposite to the applied magnetic flux, the induced magnetic flux density is calculated in a lumped way by using Ampere’s law, and the net magnetic flux density is defined by introducing the magnetic Reynolds number. In the process of the numerical integral, the skin effect is considered according to the characteristics of the eddy currents distribution. Also, the analytical results are compared with the finite element results and the causes leading to difference being analysed.


## 2. Analysis of Eddy-Current Loss

The cylindrical model of the induction device is shown in Figure2. There are n surface mounted permanent magnets radial magnetized on each rotor, with alternating pole configurations, and n is even. Where l1 and l2 are the arc length of the permanent magnets on inner rotor and outer rotor, R is the distance to the stator core center from the center of the cylindrical model; R1 and R2 are the circular radius of inner and outer permanent magnets.Figure 2
Cylindrical model.To simplify the mathematical modeling, the cylindrical model is simplified as a linear one as shown in Figure3.Figure 3
Linear model.Whereb is the axial length of model, a(a=(l1+l2)/2) is the length of permanent magnets, l is the spacing distance between adjacent permanent magnets, d is the thickness of stator, and g is the air gap. With this model simplification, the anticlockwise rotating movement at an angular velocity ω is changed into linear movement of the stator, with the constant velocity v=ωR perpendicular to the direction of alternative magnetic field B.As shown in Figure4, the movable charges in the pole projection area of the stator are subject to the magnetic force qv×B which directed upward along the stator. This action will lead to a separation of the charges, the positive ones moving to the top of the stator and the negative ones to the bottom. The Coulomb charges stack at the ends of the region of the pole projection and produce an electric field E pointing downward that will tend to decrease the total force on the given charge moving upward in the interior [9].Figure 4
Polarization in pole projection area and charge distribution in the stator core due to motion of inducing magnetic field.In the present case, the current density induced in the stator is given by(1)J=σE′,
where σ is the conductivity of the stator and the electric field intensity E′ can be calculated as [10]
(2)E′=E+Ev,
where E is the electrostatic field of Coulomb charge along the edges of the pole projection area. Ev is the electromotive field driving charges and is expressed as follows:
(3)Ev=v×B.
The current density J induced in the pole projection area can be written down by using (1), (2), and (3) as follows:
(4)J=σ(E+v×B).
Since the electromotive field v×B is easily known, E is the key step to obtain the current density J. In Figure 4,
(5)B=-Baz,v=vax.Py is defined as the polarization with respect to the y component and is written down as follows:
(6)Py=ε0Ev,
where ε0 is the permittivity in free space. Thus, Py can be expressed by using (3), (5), and (6) as follows:
(7)Py=ε0v×B=ε0vBay.
The surface charge densities on the plane parallel to the xz plane at y=±b/2 are defined as ρsy±, which is analogously calculated by using the method of calculating the magnetic pole density at the surface of the permanent magnet [11]. In the present case, ρsy+ and ρsy- are written down as follows:
(8)ρsy+=Py·ny+=ε0vBay·ny+=+ε0vB,ρsy-=Py·ny-=ε0vBay·ny-=-ε0vB,
where ny± is the unit outward normal vector on the planes parallel to the xz plane at y=±b/2. If ρly+ and ρly- are defined as the line charge densities at the points (ζ,b/2) and (ζ,-b/2), they are written down as follows:
(9)ρly+=+ε0vBΔζ,ρly-=-ε0vBΔζ,
where Δζ is the differential length in the x direction, as shown in Figure 4. If ly+ and ly- are the displacement vectors from the points (ζ,b/2) and (ζ,-b/2) to the point P(x,y), the electrostatic fields Ey+ and Ey- by the line charge densities are calculated by Coulomb’s law:
(10)Ey+=∫ρly+2πε0ly+|ly+|2=vB2π∫-a/2a/2{(x-ζ)(x-ζ)2+(y-b/2)2ax+(y-b/2)(x-ζ)2+(y-b/2)2ay}dζ,Ey-=∫ρly-2πε0ly-|ly-|2=-vB2π∫-a/2a/2{(x-ζ)(x-ζ)2+(y+b/2)2ax+(y+b/2)(x-ζ)2+(y+b/2)2ay}dζ.
Then, the total electrostatic field intensity E is obtained by using the superposition theory as
(11)E=Ey++Ey-.In order to simplify the algorithm in this case, we consider that the line charges±ε0vBΔζ at y=±b/2 extend indefinitely in the z direction and from -a/2 to +a/2 in the x direction, which is the source of E. Therefore, the z component of the electrostatic field E is ignored, and the total electrostatic field intensity E is expressed as E=Exax+Eyay, where Ex and Ey are the summations of the x and y components of Ey+ and Ey-, and they are written down as
(12)Ex=vB4πln⁡{(2x+a)2+(2y-b)2}{(2x-a)2+(2y+b)2}{(2x-a)2+(2y-b)2}{(2x+a)2+(2y+b)2},(13)Ey=vB2π[tan-12x+a2y-b-tan-12x-a2y-b-tan-12x+a2y+b+tan-12x-a2y+b].
Then, the eddy current density J induced in the stator can be obtained from (4). In the pole projection area, the x and y components of J are given as follows:
(14)Jx=σExJy=σ(Ey+vB).
By using the Lorentz force law, we have the magnetic braking force as follows:
(15)Fx=-∫vJyBdV=-d∫-a/2a/2∫-b/2b/2JyBdxdy.(16)Fy=-∫vJxBdV=-d∫-a/2a/2∫-b/2b/2JxBdxdy.In the outside of the pole projection area, the magnetic flux densityB is zero, and the eddy currents cannot produce the braking force; therefore, v is the volume of the pole projection area in (15) and (16). Since the number of poles n is even, the summation of Fy is zero. The Jy inside the pole projection is considered for the determination of Fx from (15), and we just need to calculate the Ey by (13).For the model of multipoles with alternating (NSNSNS) pole configuration, the magnetic flux of neighboring poles is symmetrical distribution of interest pole with respect toz-axis, and the surface charges at the edge of the neighboring pole projection area will affect the electric field of the interest pole. Thus, the expression of electric field intensity Ey should be revised. The electric field intensity of the interest pole can be expressed as
(17)Eyi(x,y)=vB2πK(x,y)K(x,y)=tan-12x+a2y-b-tan-12x-a2y-b-tan-12x+a2y+b+tan-12x-a2y+b.
Hence, the electric field intensity of the neighboring poles can be written down as follows:
(18)Ey(i-1)(x,y)=-Eyi(x+a+l,y)=-vB2πK(x+a+l,y)Ey(i+1)(x,y)=-Eyi(x-a-l,y)=-vB2πK(x-a-l,y).
Thus, the superposition electric field intensity inside every pole projection area is expressed as
(19)Ey=Eyi+Ey(i+1)+Ey(i-1).In (15), the original magnetic field B, produced by permanent magnets, will be decreased by the induced magnetic field generated by eddy currents. The induced magnetic flux density is calculated from the eddy current I by using Ampere’s law in order to achieve accurate results. The eddy current I is obtained from the eddy current density J by two-dimensional numerical integration. Here, we ignore the Jx, since it is much smaller than the Jy in the pole projection area, and the y component of I is expressed as
(20)Iy=2∫-a/2a/2∫0d/2Jydxdz.
Here, we ignore the skin effect in the stator. Thus, Jy is regarded as being uniform distribution on the xz plane, and the thickness of stator d is regarded as the size of the eddy-current distribution in the z-axis. Skin effect means that the current density in stator core varies exponentially inward from the surface in a time varying magnetic field. The distance in which the current density decreases to 1/e of surface value is called the nominal depth of penetration which is given by the relation
(21)δ=2ρωμ,
where ρ is the resistivity of stator core and μ is the absolute permeability of stator core. With the speed increasing, the skin effect quickly pushes the induced currents toward the surface of the stator and depletes the core of the stator from inducing eddy currents. Therefore, the thickness of the stator d is not appropriate to be regarded as the size Jy distribution in the z-axis. According to the theory of two-dimensional wave propagation in the electromagnetic field, the J in the stator is expressed as [12–14]
(22)J=Jme-αycos⁡(ωt-λx-αy)λ=πτα=σμω2=σμπnpn,
where Jm is the current density amplitude on the outer surface of stator, τ is the polar distance, σ is the conductivity of the stator, n is the number of magnetic poles, and np=60ω/2πn (r/min) is the speed. The skin effect is considered by using Jye-αz instead of Jy, and as a result the numerical integral domain of eddy current density needs to be corrected in (20). We introduce an equivalent size d′/2 instead of the constant d/2 to obtain the correct eddy current. The d′ is a function of the np and is defined as
(23)d′(np)=1-e-dαα,
where from (23), the equivalent size d′(np) decreases with the increase of speed np. Therefore, the revised Iy is a function of the np and is expressed as
(24)Iy=2∫-a/2a/2∫0d′(np)/2Jydxdz.By using Ampere’s law, the induced magnetic flux densityBin has the form of Bin=μ0Iy/2g, which decreases the original magnetic flux density B. The resultant magnetic flux density Bm is expressed as Bm=Be-Bin/B=Be-Rm, and Rm is the magnetic Reynolds number, which is defined as Bin/B. Thus, considering the induced magnetic flux, (15) is rewritten down as
(25)Fx=-∫vJyBmdV=-d∫-a/2a/2∫-b/2b/2JyBmdxdy.
Also, the eddy-current loss has the form
(26)Pe=F·v=Fxax·vax=Fxv.
Thus, we have the total eddy-current loss of the induction device given by
(27)Pe=nFxv=-ndv2BmσB×∫-a/2a/2∫-b/2b/2{12π[K(x,y)-K(x-a-l,y)-K(x+a+l,y)]+1{12π[K(x,y)-K(x-a-l,y)}dxdy.
The calculated results can be obtained through the numerical integral method. This power loss is the heat source of heat transfer fluids in the stator.


## 3. Validation by Finite Element Analysis

Two specific examples with 20 and 10 poles are given, which have the same structure parameters and physical parameters except for the number of poles. Based on these parameters listed in Table1, transient FEM models and mathematical models considering two cases of with and without skin effect are built, and the comparison between mathematical analysis results and FEM results has been carried out.Table 1
Parameters.
Description
Variable
Value
Average radius of stator (mm)
R
200
Radius of inner permanent magnet (mm)
R
1
194.5
Radius of outer permanent magnet (mm)
R
2
205.5
Radian of permanent magnets (rad)
/
π
/
60
Thickness of stator (mm)
d
5
Air gap (mm)
g
7
Axial length of model (mm)
b
200
Conductivity of the stator (1/Ωm)
σ
5.8
×
1
0
7
Absolute permeability of stator core (H/m)
μ
4
π
×
1
0
-
7
Original flux density in the air gap (T)
B
0.4The results of the two examples with 20 and 10 poles are shown in Figures5 and 6, respectively. The analytical model without considering the skin effect has small application range, and the deviation between the analytical results and the FEM results appears when the rotation speed is just over 200 r/min. However, the analytical model with considering the skin effect has wider application range, which agreed well with the FEM results within 400 r/min range of rotation speed.Figure 5
Analytical results and FEM results(n=20).Figure 6
Analytical results and FEM results(n=10).The influence of skin effect is related to the number of poles. The less the number of poles is, the smaller the frequency of AC magnetic field and the influence of skin effect are, at the same rotation speed. As shown in Figure7, the example with 10 poles has the smaller difference between two mathematical models with and without skin effect than the one with 20 poles.Figure 7
Difference between two mathematical models with and without skin effect.The analytical results have been lagged behind the finite element results, and the difference increases with the speed increase, as shown in Figure8. When the rotation speed exceeds appropriate range, the analytical model cannot provide preliminary estimations of eddy-current loss. The main reason is considered that the eddy-current loss is calculated only in the volume of the pole projection from (27), and the loss outside the pole projection is ignored which increases with the speed increase. The magnetic flux distribution at 0 r/min and 400 r/min speeds is shown in Figures 9(a) and 9(b). The available magnetic flux inside the pole projection decreases relatively with the increase of speed. On the contrary, the magnetic flux outside the pole projection increases which makes no contributions to the analytical model of eddy-current loss.Figure 8
Difference between analytical results and FEM results.Distribution of flux at speed of (a) 0 r/min and (b) 400 r/min.
(a)
(b)


## 4. Conclusion

In this paper, an induction heating device developed for heating fluid has been presented. Driven by power installation such as turbine, motor, and fan, this device can realize the conversion of the wind energy, water energy, and recovered energy of other power installations into thermal energy.A mathematical model of eddy-current loss for this device has been proposed based on the eddy-current brake concept. In modeling the eddy-current loss, the skin effect has been considered by correcting the numerical integral domain of eddy current density, which effectively extends the applied speed range of the mathematical model. Based on specific examples, this model was validated by FEM analysis. This model shows an explicit relationship between structure parameters and power performance, which is useful for design and structure optimization of the induction heating device.



---

*Source: 923745-2014-06-01.xml*