DEM Simulation of Burden Distribution in the Upper Part of COREX Shaft Furnace

Abstract

The role of COREX shaft furnace plays is similar with that of lump zone in blast furnace, but there are still some differences existed between them, especially the burden distribution. However, the research on burden distribution in COREX shaft furnace is quite few. Therefore, a three dimensional model is established in present work based on Discrete Element Method (DEM), after validated by industrial experiments, the model is used to investigate the burden profiles and distribution along radius in the upper part of COREX shaft furnace, the studied parameters include the distributor angle, rotating speed and length, and the stock line. The results show that the distributor angle affects the burden profile most, stock line is next, while rotating speed and length of distributor are least. The small particles segregate more intensively than large and medium particles. Additionally, it is better to reduce the distributor length and stock line in order to obtain a stable burden distribution along radius.

1. Introduction

COREX process is a kind of nonblast furnace smelting reducing ironmaking technology, designed for meeting resource and environmental demands.^1,2,3,4) And it has realized industrial scale in South Africa, Korea and India, especially the largest plant C-3000 at Baosteel in China.^5,6,7) According to the practical operation, it has been found that the quantity and quality indexes of COREX process are directly determined by the shaft furnace operation.^8,9,10,11) Although the role of shaft furnace plays is similar with that of lump zone in blast furnace, there are some differences existed between them in terms of gas composition and distribution, as well as burden distribution, which is of significant importance for shaft furnace.^1,12)

Many researchers including the authors pay much attention to burden distribution in blast furnace numerically and experimentally,^{13,14,15,16,17)} while few researchers studied burden distribution in COREX shaft furnace for its short history and few instruments. Li Qiang et al.¹⁸⁾ investigated the effect of distributor parameters on the radial bulk density basing on a cold model of COREX Shaft Furnace. He then studied the characteristics of binary size burdens movement based on the Discrete Element Method (DEM).¹⁹⁾ However, the model he developed is not the same size as that in practice, besides, the distributor is fixed in his model, while it is rotating during discharging burdens. Therefore, it’s of great necessity to study the burden distribution in COREX shaft furnace by establishing an actual size model with a rotating distributor.

The method adopted in present work is DEM, originated from Cundall and Strack,²⁰⁾ which has become one of the most popular and reliable simulation methods for granular behaviors. In the ironmaking field, plenty of papers about the burden distribution found that the simulation results basing on DEM coincided well with the experiments and industrial operation.^{14,15,21,22,23,24,25,26)}

In present work, a three dimensional model of charging system in COREX shaft furnace is established basing on DEM. After validated by industrial experiments, this model is used to investigate the effects of some factors on the burden profiles and size distribution along radius, such as distributor angle, rotating speed and length, and the stock line.

2. DEM Method

In DEM, it is described as inter-particle model, shown in Fig. 1, that a particle may collide with neighboring particles or wall.^20,24) The model is simplified as spring and dashpot in normal direction, spring, dashpot and slider in tangential direction. Figure 2 illustrates the interactive forces between particles. At any time τ, the governing equations for particle i can be written as follows.^25,26,27)   

$$m_i\frac{d{\mathbf{V}}_i}{d\tau }=\sum _{j=1}^K({\mathbf{F}}_{cn,ij}+{\mathbf{F}}_{dn,ij}+{\mathbf{F}}_{ct,ij}+{\mathbf{F}}_{dt,ij})+m_i\mathbf{g}$$(1)

  

$$I_i\frac{d{\omega }_i}{d\tau }=\sum _{j=1}^K({\mathbf{M}}_{t,ij}+{\mathbf{M}}_{r,ij})$$(2)

Where, m_i, V_i, and ω_i are the mass, translational velocity and angular velocity of particle i, respectively, and K is the particle number in contact with particle i. I_i is the moment of inertia of particle i, given by I_i=2/5m_iR_i². The forces involved are the gravitational force, m_ig, and inter-particle forces between particle i and j, which includes the normal and tangential contact forces, F_cn,ij and F_ct,ij, and damping forces in normal and tangential directions, F_dn,ij and F_dt,ij. Torques M_t,ij and M_r,ij arise from tangential forces (F_ct,ij + F_dt,ij) and rolling friction, respectively. All the forces and torques used in the model are listed in Table 1.^{24,25,26,28,29,30)}

Fig. 1.

The contact model in DEM.

Fig. 2.

Interactive forces between particle i and j.

Table 1. Forces and torques acting on particle i.

Forces and torques	Symbols	Equations
Normal contact force	Fcn,ij	${\mathbf{F}}_{cn,ij}=-(4/3)E_{ij}{R}_{ij}{}^{1/2}{\delta }_n{}^{3/2}\mathbf{n}$
Normal damping force	Fdn,ij	${\mathbf{F}}_{dn,ij}=-2\times {(5/6)}^{1/2}\zeta {[2m_{ij}{E}_{ij}{(R_{ij}{\delta }_n)}^{1/2}]}^{1/2}{\mathbf{V}}_{n,ij}$
Tangential contact force	Fct,ij	${\mathbf{F}}_{ct,ij}=-8G_{ij}{(R_{ij}{\delta }_t)}^{1/2}\mathbf{t}$
Tangential damping force	Fdt,ij	${\mathbf{F}}_{dt,ij}=-2\times {(5/6)}^{1/2}\zeta {[8m_{ij}{G}_{ij}{(R_{ij}{\delta }_n)}^{1/2}]}^{1/2}{\mathbf{V}}_{t,ij}$
Torque by tangential forces	Mt,ij	${\mathbf{M}}_{t,ij}={\mathbf{R}}_{ij}\times ({\mathbf{F}}_{ct,ij}+{\mathbf{F}}_{dt,ij})$
Rolling friction torque	Mr,ij	${\mathbf{M}}_{r,ij}=-{\mu }_{r,ij}|{\mathbf{F}}_{n,ij}|{\omega }_{t,ij}/\left|{\omega }_{t,ij}\right|$

Where, E_ij=[(1-υ_i²)/E_i+(1-υ_j²)/E_j]^–1, R_ij=(1/R_i+1/R_j)^–1, m_ij=(1/m_i+1/m_j)^–1, G=E/(2+2υ), G_ij=[(1-υ_i)/G_i+(1-υ_j)/G_j]^–1, ς=[lne/(π²+ln²e)^1/2]^1/2, n=δ_n,ij/|δ_n,ij|, t=δ_t,ij/|δt_,ij|, V_ij=V_j–V_i+ω_j×R_j–ω_i×R_i, V_n,ij=(V_ij·n)·n, V_t,ij=(V_ij×n)×n. Note that the tangential forces (F_ct,ij+F_dt,ij) should be replaced by the maximum friction force μ_s|F_cn,ij + F_dn,ij|t when they are larger than later force.

3. Model Simulation

Similar with the chute in blast furnace, COREX adopts GIMBAL distributor to discharge the burden from hopper to the shaft furnace circumferentially even. Since this paper mainly focuses on the burden distribution in the shaft furnace after falling from GIMBAL distributor, the model in the simulation is not the whole shaft furnace, but only the upper part with charging system. Thus the paper establishes a three-dimensional model, shown in Fig. 3, whose size is just the same as the COREX shaft furnace used in Baosteel.

Fig. 3.

Three dimensional model used in the simulation.

The height from the bottom of shaft furnace in Fig. 3 to the end of GIMBAL distributor at still condition is 5 m, which is right the stock line in practice, thus the results of simulation and industrial experiments can be easily compared under the same condition. The parameters in the simulation are the same as that of industrial operation, for example, the distributor angle between the centerlines of distributor and shaft furnace is 20°, and the rotating speed is 0.523 rad·s^–1.

In Baosteel, COREX shaft furnace is charged with mixture burden of ore, pellet, coke and flux from the hopper and distributor. As the burden sizes range widely, it will take numerous time to carry out the simulation according to the actual sizes of burden. Thus the simulation select coke as the object, for two reasons: firstly, the volumes of coke and pellet are relatively larger than other burden materials, so it’s better to choose coke or pellet; secondly, over 90% pellets’ sizes concentrate between 10 mm and 12.5 mm, it will too difficult to distinguish the size segregation clearly and the simulation results will also be unpersuasive. The detailed coke particle conditions are tabulated in Table 2, in which the sizes are adjusted relatively 2.5 times bigger than actual ones in order to reduce the calculated particles’ number.

Table 2. Coke particle conditions in the simulation.

Parameters	Large	Medium	Small	Unit
Diameter	101.88	64.38	35.00	mm
Mass ratio	27.66	53.05	19.29	wt%
Number of particles	4552	34494	78132	–

The geometries in simulated model, including hopper, GIMBAL distributor, shaft furnace and so on, are made of steel, and they are called wall generally. Table 3 presents the detailed parameters in the simulation.

Table 3. Parameters used in the simulation.

Parameters	Symbols	Particle	Wall	Unit
Density	ρ	1100	7850	kg·m–3
Shear modulus	G	2.2×106	7.9×1010	Pa
Poisson’s ratio	υ	0.22	0.30	–
Static friction coefficient	μs	0.60	0.40	–
Rolling friction coefficient	μr	0.15	0.05	–
Restitution coefficient	e	0.50	0.50	–
Time step	τ	1×10–5		s

Where, the static friction coefficients of particle and wall stand for the coefficient of particle-particle and particle-wall, respectively, so do the rolling friction coefficient and restitution coefficient.

Firstly, ternary well-mixed cokes, about 15000 kg, are generated in the hopper. Then open the gate at the bottom of hopper to discharge cokes from the hopper. And the GIMBAL distributor works at the same time. Stop the simulation at any time during discharging, like 63s, to obtain the snapshots of cross section in the GIMBAL distributor shown in Fig. 4, where L is the length of GIMBAL distributor. It can be seen that the small and medium size coke particles generally slide on the bottom of distributor, while the large ones tends to move away from the bottom. And it can be also found that the small and medium particles are likely to have greater angular velocity than the large ones for their smaller mass. These phenomena have been found by previous research,²¹⁾ therefore it can be a proof that the simulation in this paper is suitable.

Fig. 4.

Snapshots of cross section in the GIMBAL distributor at different parts (a) 0.25L, (b) 0.50L, (c) 0.75L and (d) 1.00L. (Online version in color)

In order to acquire the burden profile after discharging, a slide with certain width, through the shaft center and vertically to the horizontal plane, has been extracted. The burden profile has been measured by experiments performed in industrial furnace as well. The method is as follows, firstly, select a horizontal plane above the burden surface as the reference plane; secondly, measure the depths from the plane to the burden surface by using a plumb along radial direction in fixed interval; thirdly, convert these depths into the burden heights; at last, connect the discontinuous heights with line segments, thus the burden profile is obtained. Figure 5 gives the comparison of burden profiles between simulation and experiments, where the black line with squares and the colored particle are the results of experiment and simulation respectively. And the black, white and grey particles represent the large, medium and small cokes respectively. It can be seen from Fig. 5 that the burden profile of simulation coincides well with that of experiments. So it can be another proof of the suitable simulation.

Fig. 5.

Burden profiles of experiment and simulation.

A series of same cuboids with square bottom are extracted from the burden along plus-minus X and Y directions after discharging, and then measure the qualities of each size cokes in every cuboids. The qualities of all the cokes along the above four directions almost have no difference, thus the average qualities are selected to represent the burden qualities along radius. The comparison of coke qualities has been made between simulation and experiments, shown in Fig. 6, where the radius has been normalized.

Fig. 6.

Comparison of coke qualities between experiments and simulation.

In Fig. 6, coke qualities along radius of experiments are larger than those of simulation since the cuboids used in the experiments are bigger for containing cokes. Though the measured points of coke qualities in experiments are fewer, the tendency of experiments coincides well with that of simulation in spite of the values. As a result, the model basing on DEM has been validated and can be used to carry out further simulation.

The qualities of each size are presented in Fig. 7. The qualities of three size cokes show almost same tendencies along radius, which increase at first, and then decrease after a maximum.

Fig. 7.

Coke qualities of each size along radius.

However, the particle segregations can’t be obtained easily just according to Fig. 7. Thus the Segregation Index (SI) is induced so as to evaluate the particle size segregations clearly, and it can be expressed by Eq. (3).   

$$S{I}_k=\frac{FM{S}_k}{IM{S}_k}\times \text{100}\%$$(3)

Where, FMS and IMS represent final and initial mass fraction of each size cokes, and the subscript k means different kinds of coke, in detail, L, M and S stand for large, medium and small size cokes. If SI is bigger than 1 at a position, it means that particles tend to concentrate at this position, and vice versa.

It can be found that the SI of three kinds of coke shows different tendencies along radius from Fig. 8. Specifically, SI_L increases gradually in general along radius; SI_M rises at first, then fluctuates little in the middle, decreases at last; while SI_S is almost opposite to SI_M. This is because that large and medium particles locate further from the from the centerline than small ones during slide from the distributor, which results in the same phenomenon in the burden stream, thus these particles are relatively easy to roll toward the wall. It can be also seen that SI_S fluctuates heavier than SI_L and SI_M, so the standard deviation is adopted to calculate the segregation degree. The standard deviations from large to small particles are 0.162, 0.142 and 0.520, which implies that small particles segregate heaviest, and large and medium particles segregate relatively lighter. The reason is that small particles strike to the burden close to the center, and also they penetrate easily between large and medium particles.

Fig. 8.

SI of each size coke along radius.

4. Results and Discussion

In consideration of practical operation and some adjustable parameters in actual COREX shaft furnace, the model is used to study the influence of some factors, like distributor angle, rotating speed and length, and stock line, on burden profiles and particle size segregations. Table 4 lists the factors and values investigated in the simulation.

Table 4. Factors and values investigated in simulation.

Factors	Values	Unit
Distributor angle	10°, 20° and 30°	°
Rotating speed	0.523, 1.047 and 1.570	rad·s–1
Distributor length	L, 2L and 3L	m
Stoke line	4, 5 and 6	m

Where, the length of GIMBAL distributor is short for L due to the confidentiality agreements with Baosteel.

4.1. Distributor Angle

Figure 9 shows the burden profiles at different distributor angles. When the distributor angle changes from 10° to 30°, the position of burden apex moves towards the furnace wall, and the height of burden apex decreases. This is because that the main burden streams strike further, specifically, the striking points are about 0.32 R, 0.50 R and 0.68 R, respectively, thus the burden apex heights reduce for the fixed burden volume.

Fig. 9.

Burden profiles at the distributor angles of (a) 10°, (b) 20° and (c) 30°.

Figure 10 illustrates the burden qualities and SI along radius. Since burden profile and coke qualities along radius coincide well with each other, the coke quality line along radius in Fig. 10 is to provide a reference for the explanation of SI.

Fig. 10.

Coke qualities and SI along radius at the distributor angles of (a) 10°, (b) 20° and (c) 30°.

In Fig. 10(a), the SI_M tendency along the radius stays changeless but a little decrease, however, the SI_L and SI_S tendencies change greater, in detail, the SI_S along radius is almost “U”-shape like, while the SI_L tendency shows an opposite tendency to SI_S. The SI_L bigger than 1.0 is left to the apex of coke quality, which means that the large cokes tend to roll towards furnace wall. On one hand, SI_S is large near the furnace center because the small cokes are at the inner side of burden stream and they can percolate between other particles. On the other hand, it is large near the furnace wall for the hit of large and medium cokes, the percolation of small particles as well. The standard deviations of SI_L, SI_M and SI_S are 0.339, 0.104 and 0.351, respectively, which also proves that the segregations of large and small cokes are heavier than medium ones.

In Fig. 10(c), the SI_L and SI_M show similar tendency along radius, which go up initially and then stay unchanged, and decrease a little finally. This means that the large and medium cokes concentrate around the apex of burden quality line. The phenomenon is caused by that the room is not enough for them to roll between burden stream and wall. While SI_S generally reduces at first and then stays invariant, and rises a little in the end. This is because the small cokes are compelled to percolate and roll towards the inner part for this reason. As for the first part of SI_S tendency, it rise a little to the largest, this can be explained that the total coke qualities there are so small that few particles affect great for SI. In the terms of standard deviation, they are 0.331, 0.123 and 0.699 from large to small cokes. It can be concluded that the segregations of small cokes are the heaviest. The SI standard deviations at different distributor angles are tabulated in Table 5.

Table 5. Standard deviations of SI at different distributor angles.

Distributor angle	SIL	SIM	SIS
10°	0.339	0.104	0.351
20°	0.162	0.142	0.520
30°	0.331	0.123	0.699

It can be found that the standard deviation of SI_S is the largest while SI_M is the smallest no matter how many the angles are. When the distributor angle increases, the standard deviations of SI_L go down and then rise, the SI_M changes little, while the SI_S rises greatly. So it can be confirmed that the segregations of small cokes grow severe with the angle rising, and the segregations of large cokes are relatively slight when the angle is not too large or small.

4.2. Distributor Rotating Speed

The burden profiles, shown in Fig. 11, are nearly the same at different rotating speeds. So the coke qualities along radius are induced in order to investigate whether there exits some differences or not, and the comparison result is given in Fig. 12.

Fig. 11.

Burden profiles at the rotating speeds of (a) 0.523 rad·s^–1, (b) 1.047 rad·s^–1 and (c) 1.570 rad·s^–1.

Fig. 12.

Coke qualities along radius at different rotating speeds.

It can be found that the apexes of coke quality line transfer to right a little along the radius, and the coke qualities values of the apex decline from 51.4 kg to about 47.2 kg with the distributor rotating speed increasing. Meanwhile, the coke qualities near the wall side are from 11.4 kg to 14.1 kg, showing an increase tendency. All the above phenomena imply that the burdens are likely to strike a little further while the distributor rotating speed rises.

Figure 13 provides the burden qualities and SI along the radius. The SI tendencies are almost the same except for few points comparing with the base model whose rotating speed is 0.523 rad·s^–1. So it suggests that the distributor rotating speed has little effect on coke particle size segregations.

Fig. 13.

Coke qualities and SI along radius at the rotating speeds of (a) 0.523 rad·s^–1, (b) 1.047 rad·s^–1 and (c) 1.570 rad·s^–1.

Table 6 lists the standard deviations of SI at different rotating speed. The standard deviations of SI don’t change much with the rotating speed varying. It shouldn’t be ignored that the standard deviation of SI_S is the largest, which implies that small cokes segregate much in spite of rotating speed.

Table 6. Standard deviations of SI at different rotating speeds.

Rotating speed	SIL	SIM	SIS
0.523 rad·s–1	0.188	0.113	0.465
1.047 rad·s–1	0.162	0.142	0.520
1.570 rad·s–1	0.170	0.136	0.510

4.3. Distributor Length

The burden profiles at different distributor lengths are presented in Fig. 14. There seems to be no difference among the three burden profiles. Thus the burden quality lines along radius are used to check it quantitatively in Fig. 15.

Fig. 14.

Burden profiles at the distributor lengths of (a) L, (b) 2L and (c) 3L.

Fig. 15.

Coke qualities along radius at different distributor lengths.

Generally, the burden quality profiles shift towards wall direction slightly, and the left half moves greater than the right half. This is because that the height between distributor and burden surface becomes shorter when the distributor length rises, which leads the burden to strike further from the center. While the velocities become smaller, which means that the energies are not great enough to roll further for the right half burden.

Figure 16 shows the burden qualities and SI along radius. The tendencies of SI at different distributor lengths are nearly the same, specifically, SI_L and SI_M increase at first 0.3 R and then fluctuate about 1 along radius, SI_S tendency seems like a “U” shape. The calculated standard deviations of SI are list in Table 7.

Fig. 16.

Coke qualities and SI along radius at the distributor lengths of (a) L, (b) 2L and (c) 3L.

Table 7. Standard deviations of SI at different distributor lengths.

Distributor lengths	SIL	SIM	SIS
L	0.162	0.142	0.520
2L	0.194	0.141	0.534
3L	0.258	0.146	0.636

It can be found that the standard deviations of SI_L and SI_S are getting larger while the SI_M almost keeps unchanged with the increase of distributor length. This is the result of narrower burden stream width, which forms severe size segregations along the whole radius. Hence it implies that shorter distributor length can help to reduce the segregations.

4.4. Stock Line

In terms of stock line, it can be simply regarded as the height between the end of distributor and the burden surface. Figure 17 displays the burden profiles at different stock lines.

Fig. 17.

Burden profiles at the stock lines of (a) 4 m, (b) 5 m and (c) 6 m.

When the stock line goes up, the burden profile apex is likely to move towards the wall direction, and the height declines at the same time. In detail, the positions of burden apex are about 0.45 R, 0.50 R and 0.59 R, respectively. The main reason is that the burden streams strike further from the furnace center since the burden discharged from the distributor travels longer with the stock line rising.

Figure 18 shows the burden qualities and SI along radius at different stock lines. For all the cases of stock line, the SI_L displays the similar tendency along radius, so do the SI_M and SI_S. So it infers that the stock line has little influence on the SI tendency. However, the segregation degree is not the same quantitatively, the results of standard deviations are given in Table 8.

Fig. 18.

Coke qualities and SI along radius at the stock lines of (a) 4 m, (b) 5 m and (c) 6 m.

Table 8. Standard deviations of SI at different stock lines.

Stock lines	SIL	SIM	SIS
4 m	0.147	0.118	0.440
5 m	0.162	0.142	0.520
6 m	0.280	0.109	0.582

It can be found that SI_L and SI_S increase with the rise of stock line, which means that the segregations of large and small cokes become more intensive. So it’s better to decline the stock line in order to obtain the relatively uniform distribution along radius.

5. Conclusions

The paper establishes a three dimensional model validated by experiments in order to study the burden profiles and coke size segregations along radius in the upper part of COREX^® shaft furnace basing on DEM. The factors simulated in the paper include GIMBLE distributor angle, rotating speed and length, and stock line. Generally, the smaller the particle is, the larger angular velocity it has. And the small particles are likely to press to the GIMBLE distributor bottom. These phenomena mean that the segregations of angular velocity and particle size happen on the GIMBLE distributor, which directly results in the burden segregations in shaft furnace. Additionally, the small particles segregate more intensively than large and medium particles in shaft furnace.

In terms of the ranges of factors investigated in the paper, the distributor angle affects the burden profiles most greatly, the next is stock line, while the rotating speed and length of distributor affect least. When the factors increase, the standard deviations of SI_S increase for distributor angle, length and stock line, almost stay unchanged for distributor rotating speed; while the standard deviations of SI_L firstly decrease and then increase for distributor angle, fluctuate a little for distributor rotating speed, increase for distributor length and stock line; the standard deviations of SI_M fluctuate a little for all the factors. Therefore, in order to achieve a uniform burden distribution along radius and smooth operation in COREX shaft furnace, it is recommended to decline the length and stock line, keep the distributor angle at 20° position.

This work is an initial explore on burden distribution in COREX shaft furnace by using DEM, further simulation will be carried out for the consideration of the irregular particle shape. And the simulation will also take the effect of gas flow, particle softening and shrinkage, multiphase flow including thermo-chemical behavior into account inside the entire shaft furnace or even the COREX process. Work in this direction is undertaken and will be reported hopefully in the future.

Nomenclature

e: restitution coefficient (–)

E: Young’s modulus (Pa)

F: force between particles (N)

g: gravitational acceleration (m·s^–2)

G: shear modulus (Pa)

I: moment of inertia (kg·m^–2)

K: particle number in contact with particle i

m: mass of particle (kg)

M: torque on particle (N·m)

n: unit vector in normal direction

R: particle radius (m)

t: unit vector in tangential direction

V: velocity (m·s^–1)

Greek letters

δ: displacement (m)

ζ: damping ratio (–)

μ: friction coefficient (–)

ρ: density (kg·m^–3)

τ: time (s)

υ: Poisson’s ratio (–)

ω: angular velocity (rad·s^–1)

Subscripts

c: contact

d: damping

i: ith particle

ij: between particle i and j

j: jth particle

n: normal component

r: rolling

s: static

t: tangential component

k: kind of coke size

Acknowledgement

The authors would like to appreciate much for the anonymous reviewers and editors for the improvement of this work, and Mr. Jian Xu from Chongqing University for suggestions on this work and Miss Daichun Wei for language review.

References

• 1)   S.  Wu,  J.  Xu,  S.  Yang,  Q.  Zhou and  L.  Zhang: ISIJ Int., 50 (2010), 1032.
• 2)   S.  Wu,  Q.  Zhou,  J.  Xu,  S.  Yang,  Y.  Du and  L.  Zhang: Proc. of 5th Int. Cong. on the Sci. Technol. of Ironmaking, Chinese Society for Metal, Beijing, (2009), 1280.
• 3)   X.  Liu,  G.  Pan,  G.  Wang and  Z.  Wen: Energ. Fuel., 25 (2011), 5729.
• 4)   E.  Ottenschlaeger,  W.  Kepplinger,  G.  Papst and  R.  Hauk: MPT, 9 (1986), 24.
• 5)  Steel Times Int., 13 (1989), 24.
• 6)   A.  Eberle,  E.  Eichberger,  C.  Bohm,  G.  Leonard and  H.  Wetzling: Iron Steel Eng., 75 (1998), 25.
• 7)   A.  Eberle,  D.  Siuka and  C.  Böhm: Stahl Eisen, 126 (2006), 31.
• 8)   H.  Xu,  L.  Guo,  Z.  Zhang and  J.  Zhu: Baosteel Tech. Res., 5 (2011), 44.
• 9)   W.  Xu,  Y.  Guo and  C.  Wang: Baosteel Tech. Res., 5 (2011), 45.
• 10)   K.  Luo and  M.  Jin: Ironmaking, 30 (2011), 62.
• 11)   Q.  Zhang,  L.  Guo and  X.  Chen: Proc. of 5th Int. Cong. on the Sci. Technol. of Ironmaking, Chinese Society for Metal, Beijing, (2009), 1230.
• 12)   W.  Ying,  Y.  Sun,  Z.  Luo and  Z.  Zou: China Metall., 22 (2012), 13.
• 13)   J.  Xu,  S.  Wu,  M.  Kou,  L.  Zhang and  X.  Yu: Appl. Math. Model., 35 (2011), 1439.
• 14)   H.  Mio,  S.  Komatsuki,  M.  Akashi,  A.  Shimosaka,  Y.  Shirakawa,  J.  Hidaka,  M.  Kadowaki,  S.  Matsuzaki and  K.  Kunitomo: ISIJ Int., 48 (2008), 1696.
• 15)   H.  Mio,  S.  Komatsuki,  M.  Akashi,  A.  Shimosaka,  Y.  Shirakawa,  J.  Hidaka,  M.  Kadowaki,  H.  Yokoyama,  S.  Matsuzaki and  K.  Kunitomo: ISIJ Int., 50 (2010), 1000.
• 16)   T.  Inada,  Y.  Kajiwara and  T.  Tanaka: ISIJ Int., 29 (1989), 761.
• 17)   M.  Kou,  S.  Wu,  J.  Xu,  X.  Guo,  C.  Fu and  K.  Du: SCANMET IV - 4th Int. Conf. on Proc. Dev. in Iron and Steelmaking, MEFOS, Lulea, (2012), 441.
• 18)   Q.  Li,  D.  Liu,  L.  Zhang and  Z.  Zou: J. Northeast. Univ. (Nat. Sci.), 33 (2012), 542.
• 19)   Q.  Li,  D.  Liu,  L.  Zhang and  Z.  Zou: J. Northeast. Univ. (Nat. Sci.), 33 (2012), 848.
• 20)   P. A.  Cundall and  O. D. L.  Strack: Geotechnique, 29 (1979), 47.
• 21)   H.  Mio,  S.  Komatsuki,  M.  Akashi,  A.  Shimosaka,  Y.  Shirakawa,  J.  Hidaka,  M.  Kadowaki,  S.  Matsuzaki and  K.  Kunitomo: ISIJ Int., 49 (2009), 479.
• 22)   H.  Mio,  M.  Kadowaki,  S.  Matsuzaki and  K.  Kunitomo: Miner. Eng., 33 (2012), 27.
• 23)   S.  Natsui,  S.  Ueda,  Z.  Fan,  N.  Andersson,  J.  Kano,  R.  Inoue and  T.  Ariyama: ISIJ Int., 50 (2010), 207.
• 24)   J.  Zhang,  Z.  Hu,  W.  Ge,  Y.  Zhang,  T.  Li and  J.  Li: Ind. Eng. Chem. Res., 43 (2004), 5521.
• 25)   Y.  Yu and  H.  Saxén: Chem. Eng. Sci., 65 (2010), 5237.
• 26)   Y.  Yu and  H.  Saxén: ISIJ Int., 52 (2012), 788.
• 27)   Z. Y.  Zhou,  H. P.  Zhu,  A. B.  Yu,  B.  Wright and  P.  Zulli: Comput. Chem. Eng., 32 (2008), 1760.
• 28)   A.  Di Renzo and  F. P.  Di Maio: Chem. Eng. Sci., 59 (2004), 525.
• 29)   Y. C.  Zhou,  B. D.  Wright,  R. Y.  Yang,  B. H.  Xu and  A. B.  Yu: Physica A, 269 (1999), 536.
• 30)   Y.  Tsuji,  T.  Tanaka and  T.  Ishida: Powder Technol., 71 (1992), 239.