One-dimensional model of heat-recovery, non-recovery coke ovens

JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 Fuel xxx (2016) xxx–xxx 1 Contents lists available at ScienceDirect Fuel journal homepage: www.elsevier.com/locate/fuel 5 6 4 One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model 7 Rafal Buczynski a,⇑, Roman Weber a, Ronald Kim b, Patrick Schwöppe b 3 8 9 a b Institute of Process Energy Engineering and Fuel Technology, Clausthal University of Technology, Agricola Str. 4, 38678 Clausthal-Zellerfeld, Germany ThyssenKrupp Industrial Solutions AG, Business Unit Process Technologies, Friedrich-Uhde-Str. 15, 44141 Dortmund, Germany 10 11 h i g h l i g h t s 1 3 14 15  Moving boundary technique is used to track moisture evaporation and condensation in coking-beds. 16  Carbonization model provides both yield and composition of raw gas at any instance of the process. 17  Variations of the fixed-bed properties (including coking pressure) during carbonization have been modelled. 18 a r t i c l e 2 1 0 3 21 22 23 24 25 26 27 28 29 30 i n f o Article history: Received 18 November 2015 Received in revised form 25 January 2016 Accepted 27 January 2016 Available online xxxx Keywords: Heat-recovery coke ovens Carbonization Moving boundary method a b s t r a c t A one-dimensional mathematical model of HR/NR coke ovens has been developed and it includes a series of sub-models. The heart of the model is a hydraulic network sub-model described in Part I (Buczynski et al., submitted for publicaion). In this paper (Part II) the carbonization process is described and casted into a coking-bed sub-model. The sub-model handles heat transfer, moisture evaporation and condensation as well as devolatilization as time-dependent processes progressing along the bed height. The fixed bed properties like the effective thermal conductivity, porosity and bulk density vary with time and degree of carbonization. Novelty of the work is in application of the moving boundary technique in solving the heat balance equation in order to remedy the discontinuities occurring due to moisture evaporation and condensation. Ó 2016 Published by Elsevier Ltd. 32 33 34 35 36 37 38 39 40 41 42 43 44 45 1. Introduction and objectives 46 As it has been formulated in our previous publication [1], the overall objective of our work is the development of a onedimensional time-dependent mathematical model of heatrecovery non-recovery coke ovens. The model is to serve as a tool for both analysis and optimization of horizontal ovens. While our Part-I publication [1] provides both a general description of the model and a detailed description of the hydraulic network submodel, the current paper is concerned with mathematical description of the processes proceeding in the coal/coke fixed bed. The fixed bed model is named here as the coking-bed sub-model. There exists a wealth of literature concerning mathematical description of processes progressing in fixed-beds and carbonization is one example only. Studies dealing with fixedbed combustion of coal appeared already in the 1950’s [2]. More sophisticated descriptions of both combustion and gasification of 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ⇑ Corresponding author. coals in fixed-beds can be found for example in [3–9]. More recent publications have appeared which are related to mathematical modelling of fixed-bed combustion for various types of biomass [9–14] and refuse derived fuels [15]. Fixed-bed sub-models are parts of larger software packages used for performance optimization of small scale domestic stoves [16–21] as well as industrial grate stokers [22–24]. Since in our Part I paper [1] a comprehensive review concerning modelling of coke ovens have been presented, in this publication we cite works directly relevant to the current paper only; these are the publications of Merrick et al. [25–29] and Klose and Nowack [30]. It will be shown later that our method of modelling of both moisture evaporation and condensation is substantially different from the Merrick’s approach [28] whilst coal pyrolysis is modelled as proposed by Merrick [25]. Klose and Nowack [30] present a two-dimensional model of coking process proceeding in a horizontal chamber reactor. The model is used to predict fields of coal temperature and flow rate of gas released during fuel decomposition. The sub-model describing heat conduction in the porous media consider phenomena such as: solid heat conduction, gas E-mail address: rafal.buczynski@tu-clausthal.de (R. Buczynski). http://dx.doi.org/10.1016/j.fuel.2016.01.086 0016-2361/Ó 2016 Published by Elsevier Ltd. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 2 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Nomenclature Greek b Dmi Dt dx Dx Dxevap e e ¼ e1 j qbulk qtrue s shape parameter mass of the daf coal (kgraw/kgdaf) time step (s) distance between grid points (m) control volume size (m) distance travelled by evaporation front (m) coal/coke bed porosity (m3gas /m3bed Þ e2 width parameter frequency factor (1/s) bulk coal/coke density (kg/m3) true coal/coke density (kg/m3) time (s) Variables A A a b c cbulk cbulk ck;bulk E E0 F0 (E) f f 1 gi gm h k keff,bulk m moðsÞ _ cond m _ Bcond m _ Tcond m mdev _ vol m _ Bdev m _ Tdev m _ evap m mi mi(E) i m ðsÞ mi Mm n n o p pcoking pboil sat pTsat q cross-sectional area (m2) matrix of constants mass fraction of ash (kgash/kgbulk) vector of constants mass fractions of carbon in the daf coal (kg/kgdaf) coal/coke charge specific heat (J/kg K) mean specific heat of the entire coal bed (J/kg K), mean specific heat of the coal/coke located in the k-th numerical cell (J/kg K) energy activation (J/kmol) ‘‘starting” activation energy (J/kmol) non-cumulative Rosin–Rammler function previous time step current time step mass fraction of the i-th raw-gas component (kgi/kgraw) mass fraction of moisture in the coking-bed (kgm/kg) mass fractions of hydrogen in the daf coal (kg/kgdaf) thermal conductivity (W/mK) coal/coke effective thermal conductivity (W/mK) vector of final yields of coke and volatiles (kgi/kgdaf) mass of char remaining in the charge at time s (kgchar) mass flow of condensating moisture (kgm/s) mass flow of condensating moisture (moving upwards) (kgm/s) mass flow of condensating moisture (moving downwards) (kgm/s) amount of volatiles released (kg) mass flow rate of raw gas (kgvol/s) mass flow of volatile matter (released at the bottom) (kgvol/s) mass flow of volatile matter (released at the top) (kgvol/s) mass flow rate of evaporating moisture (kgm/s) final yields of coke and volatiles (kgi/kgdaf) cumulative volatile release at time t (kgi/kgdaf) final yields of coke and volatiles (kgi/kgdaf) cumulative volatile release at time s (kgi/kgdaf) molar mass of moisture (kgm/kmol) mass fractions of nitrogen in the daf coal (kg/kgdaf) number of grid elements mass fractions of oxygen in the daf coal (kg/kgdaf) ASTM volatile matter content (kg/kgdaf) mean coking pressure (kPa) saturated vapour pressure at boiling temperature (N/m2) saturated vapour pressure at cell temperature (N/m2) heat flux (W/m2) universal gas constant (J/kmol K) enthalpy of evaporation (J/kgm) mass fractions of sulphur in the daf coal (kgi/kgdaf) linearized source term (W/m3) heat sources due to steam condensation (W/m3) heat sources/sinks due to thermal decomposition of coal (W/m3) heat sinks due to moisture evaporation (W/m3) Sevap explicit term in the source term linearization Sex Sim implicit term in the source term linearization T temperature (K) boiling temperature (K) Tboil coal/coke charge temperature (K) Tbulk T o ¼ 298:15 K the reference temperature (K) T k;bulk coal/coke temperature in the k-th cell (K) w mass fraction of moisture (kgm/kgbulk) V cell volume (m3) V mass fractions of volatile matter in the daf coal (kgi/ kgdaf) w weights wevap evaporation front velocity (m/s) wBevap evaporation front velocity (the bottom front) (m/s) wTevap evaporation front velocity (the top front) (m/s) xevap position of evaporation front (m) xTevap position of evaporation front moving downwards (m) position of evaporation front moving upwards (m) xBevap mass fractions of ash in the deposit (kga/kgcoal) xa xcp mass fraction of combustible part in the deposit (kgcp/kgcoal) mass fractions of moisture in the deposit (kgm/kgcoal) xm mass fraction of i-th char component (kgi/kgchar) yi R r evap s S Scond Sdev Subscripts a ash atm atmospheric boil boiling bulk charge coking coking cond condensation cp combustible part (coke) d dry coal dev raw pyrolysis gas eff effective evap evaporation ex explicit im implicit m wet coal m moisture sat saturation true true vol volatile matter (raw gas) Superscripts B bottom sub-domain boil boiling point f current time step (new, unknown value ) f 1 previous time step (old, known value) j volatile matter component T top sub-domain T cell temperature Top boundary condition at the charge top ðsÞ current time step ðs þ DsÞ next time step Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 3 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Rest B BC daf N, P, S 81 82 83 84 85 86 87 heat conduction, radiation in the fissure volume of the semi coke filled with absorbing gas, and convective transport by flowing gas. Distribution of gas in the coking charge is estimated using a modified permeability tensor of the anisotropic matrix. It is assumed that the flow resistance depends on gas velocity and flow direction. In addition, the evaporation and condensation of water and tar have been modelled. 88 2. Coking-bed sub-model 89 Fig. 1 shows a simplify sketch of the coal/coke bed. From the top, the bed is heated by radiation and convection while at the bottom heat is conducted through the sole-flue ceiling. Coke making is inherently time dependent process with heat transfer through the bed determining the temperature distribution along the bed height. In hot regions of the bed, moisture evaporates and condensates in cold regions. Pyrolysis of the coal takes place in the regions where the temperature is high enough and pyrolysis gases (named here as raw-gas) are released to the space above the bed (upperoven). Both the amount and composition of the pyrolysis gases are a function of coking time. Heat conduction through the coking-bed is calculated as proposed by Merrick [28]: The heat sink (see Eq. (1)) due to evaporation of moisture at boiling temperature Tboil is calculated by: Sevap 90 92 93 94 95 96 97 98 99 100 101 102 104 qbulk cbulk   @T bulk @ @T bulk keff;bulk ¼ þ Sevap þ Scond þ Sdev @x @s @x ð1Þ 107 where qbulk – bulk coal/coke density (kg/m3), cbulk – coal/coke charge specific heat (J/kg K), Tbulk – coal/coke charge temperature (K), keff,bulk – coal/coke effective thermal conductivity (W/mK), Scond 108 – the heat sources due to steam condensation (W/m3), Sevap – the 105 106 109 110 heat sinks due to moisture evaporation (W/m3), Sdev – the heat sources/sinks due to thermal decomposition of coal (W/m3). _ evap r evap
m ¼ V ð2Þ _ evap – mass flow rate where r evap – enthalpy of evaporation (J/kgm), m of evaporating moisture (kgm/s), V – cell volume (m3). The heat source (see Eq. (1)) due to condensation of water vapour in the cold areas of the charge is determined by: Scond 91 control volume faces (North, South) the top part of the deposit (crown side) the moisture evaporated in the top sub-domain and condensed in the central sub-domain n, s T TC the bottom part of the deposit (sole-flue side) moisture evaporated in the bottom sub-domain dry ash free grid points (North, Centre, South) _ cond r evap
m ¼ V ð3Þ _ cond – mass flow of where r evap – enthalpy of evaporation (J/kgm), m condensating moisture (kgm/s), V – cell volume (m3 ). The additional heat source/sink Sdev (see Eq. (1)) takes into account the thermal effects associated with thermal decomposition of the coal blend. These effects are described in Section 3.6. In order to solve Eq. (1) the coking-bed is divided into I elements of 15 m  4 m  dxvolume, as shown in Fig. 1. Eq. (1) is discretized in space (x-coordinate) and time. For the space discretization a central differencing is used while a fully implicit scheme is applied for discretization in time. Discretization of Eq. (1) [31] leads to:  2  f  kn T N T Pf f 1 TP ¼4 ðdxÞn Dx  f ðqcÞP T Dt P þ Sim DxT Pf  ks T Pf ðdxÞs 111 112 113 115 116 117 118 119 120 122 123 124 125 126 127 128 129 130 131 132 133 134 3 f TS 5 þ Sex Dx ð4Þ 136 where q ¼ qbulk – density (kg/m3), c ¼ cbulk – specific heat (J/kg K), k ¼ keff;bulk – thermal conductivity (W/mK), T = Tbulk – temperature (K), N, P, S – indicate grid points shown in Fig. 2, n, s – control volume faces, Sex, Sim – quantities from source term linearization, dx – distance between grid points (m), Dx – control volume size (m), Dt – time step (s), f – indicates current time step (new, unknown value ), f 1 – previous time step (old, known value). The resulting working matrix equations are: 137 aP T Pf 147 ¼ aN T Nf þ aS T Sf þb ð5Þ where 138 139 140 141 142 143 144 145 148 149 kn aN ¼ ðdxÞn ð6Þ 151 152 ks aS ¼ ðdxÞs Fig. 1. Coking-bed and its one-dimensional numerical representation. ð7Þ Fig. 2. Internal (left) and boundary (right) control volumes. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 154 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 4 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 155 157 1 b ¼ Sex Dx þ a0P T fP ð8Þ where ð9Þ ks aS ¼ 2 ðdxÞs 158 160 161 163 164 165 166 168 169 170 171 173 174 175 176 178 a0P ¼ ðqcÞP Dx Dt 182 183 184 186 187 188 190 Sim Dx ð10Þ The thermal conductivity for interface ‘‘n” is calculated as the harmonic mean of N and P values:    2kP T Pf kN T Nf   kn ¼   kP T Pf þ kN T Nf ð11Þ The thermal conductivity for interface ‘‘s” is calculated in the same way but using grid points P and S:     2kS T Sf kP T Pf   ks ¼   kS T Sf þ kP T Pf ð12Þ The density and the specific heat at interfaces ‘‘n” and ‘‘s” are determined by the following equations: ðqcÞPs ¼ 1 T Pf Z T Pf Z T Nf 1 T Sf f TN f qðT ÞcðT ÞdT ð13Þ TP f TS T Pf qðT ÞcðT ÞdT ð14Þ The product of density and specific heat at grid point P is calculated as: ðqcÞP ¼ wPs ðqcÞPs þ wPn ðqcÞPn ð15Þ where wPs and wPn are weights defined below: wPn ¼ 1 ðdxÞn 1 ðdxÞn þ ðdx1Þ ð16Þ s 191 193 194 195 197 wPs ¼ 1 ðdxÞs 1 ðdxÞn þ ðdx1Þ ð17Þ s The source terms (see Eq. (1)) have been linearized as follows: 1 S ¼ Sf þ 198  @S @T f 1 T Pf T fP 1  ð18Þ 200 S ¼ Sex þ Sim T Pf 201 where the explicit and the implicit parts are: 202 204 Sex ¼ S  @S @T f 1 205 207 208 209 210 211 Sim ¼  @S @T f 214 215 217 ð19Þ f 1 T fP 1 ð20Þ b ¼ Sex Dx þ a0Top T fTop1 2qTop ð25Þ 221 ð21Þ When the boundary condition of the second kind (given heat flux) is used at the top boundary of the computational domain, the following discretization equation is used (see Fig. 2, right):  2 T fTop1 ¼ 24 qTop  f ks T Top ðdxÞs f þ Sex Dx þ Sim DxT Top 3 f TS 5 ð22Þ In the matrix system the above equation can be written as f aTop T Top ¼ aS T Sf þ b 224 225 a0Top ðqcÞTop Dx ¼ Dt ð26Þ aTop ¼ aS þ a0Top Sim Dx 227 ð27Þ 230 Implementation of the boundary condition at the bottom boundary of the computational domain is carried out in the same way. The system of I-1 equations is then solved at any instance with the time step Dt using the tridiagonal matrix algorithm (TDMA). When at the domain boundary, the temperature is specified (boundary conditions of the first kind) no additional equations are required and the system of equations is solved to determine I-2 temperatures using the TDMA algorithm. As the coking time elapses, moisture evaporation fronts move from the top (see Fig. 3) and from the bottom towards the bed centre. This requires a special front tracking technique described in the subsequent section. 231 2.1. Moisture evaporation 243 Both moisture evaporation and condensation exert a strong influence on devolatilization rate and on the entire coke making process because distribution/flow of the moisture in the coal/coke bed affects the charge (fixed-bed) heating rate. Evaporation front moves slowly from the bottom (the sole-flue side) and from the top (the crown side) toward the centre of the deposit. Moisture evaporation precedes devolatilization. Plastic phase, created nearby devolatilization zones, forms a partially or fully impermeably barrier for the evaporated moisture so that it cannot leave the bed. Instead the moisture diffuses through the voids towards the bed centre and condenses in the cold regions. As soon as the two condensation fronts merge together (somewhere in the centre of the bed), the moisture cannot condensate anymore and leaves the coal bed through the pores in the bed or near the walls of the coke oven. There are two evaporation fronts: one at the top ( xTevap – moving 244 downwards) and the second one at the bottom ( xBevap – moving upwards); Fig. 3 shows the top front moving with the velocity wTevap . Velocities of the internal boundaries (evaporation fronts) are calculated using energy balances written for the top front: 260 @T d d keff;bulk bulk @xd @T m m keff;bulk bulk @xm ¼ qbulk
g m
r evap
wTevap ð28Þ and for the bottom one 1 Dx  f ðqcÞTop T Dt Top 213 ð24Þ 228  179 181 219 222 a0P aP ¼ aN þ aS þ ðqcÞPn ¼ 218 ð23Þ @T m m keff;bulk bulk @xm @T d d keff;bulk bulk @xd 233 234 235 236 237 238 239 240 241 242 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 261 262 263 264 265 267 268 269 ¼ qbulk
g m
r evap
wBevap ð29Þ where keff,bulk – effective thermal conductivity of the coal/coke bed (W/mK), revap(Tboil) – enthalpy of evaporation (kJ/kgm), qtrue – true coal/coke density (kg/m3), gm – mass fraction of moisture in the coking-bed (kgm/kg), wTevap and wBevap – evaporation fronts velocities (m/s) for the top and the bottom fronts, respectively. Thus, it is assumed that the evaporation rate is directly proportional to the rate of coal heating. The distance travelled by the evaporation front may be calculated using: Dxevap ¼ wevap Ds 232 ð30Þ Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 271 272 273 274 275 276 277 278 279 280 282 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 5 Fig. 3. Top evaporation front moving downwards with wTevap velocity. Fig. 4. Snapshot of the computational (coking) domain with division into top, central and bottom sub-domains. 283 284 285 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 where Ds is the time step (s) while the mass flow rate of evaporated moisture is calculated as: _ evap ¼ A
wevap
qtrue
ð1 m Þ
g m ð31Þ where A is the cross sectional surface area (m2) and e is the coal/coke bed porosity (m3gas /m3bed Þ. To accurately calculate the heating rate of the bed it is necessary to use a special numerical method to track the evaporation fronts. The discontinuity, which occurs at the evaporation fronts, is a non-trivial problem for numerical methods. To remedy the discontinuity, the computational domain is divided (along the cokingbed height) into three connected sub-domains which are named as the top sub-domain, the central sub-domain and the bottom subdomain, as shown in Fig. 4. The split of the fixed bed into the three sub-domains is based on the boiling temperature of water (Tboiling); in both the top and bottom sub-domains the (charge) temperatures are larger than Tboiling so there is no moisture present since it has already being given off. Consequently, if the temperature of these sub-domains is high enough, the pyrolysis may occur. In the central sub-domain the moisture evaporates at Tboiling and it condensates when the temperature is lower so that a temperature profile as shown in Fig. 4 is applicable. Obviously, at the beginning of the process (at t = 0), both the top and bottom sub-domains are infinitely thin (or, in other words, they do not exist) and their thickness increases with time until the evaporation fronts meet. The top sub-domain is linked with the central one using the moving boundary condition of the 1st kind (prescribed temperature). In the same way the central sub-domain is linked with the bottom one. Thus, the moving boundary temperature is always equal to the boiling temperature (for 1 bar equals 373 K, for 1.5 bar equals 385 K). In Section 2, we have described the numerical scheme for solving Eq. (1) using either the boundary condition of the second kind (heat flux given) or the first kind (temperature given). Now, we describe how the moving boundary technique is used to track the evaporation fronts. At the beginning of coking process (at time zero), the initial coal bed temperature is uniform throughout the bed. We divide the bed into I-1 numerical elements (slices) and prescribe, say one third of the elements to the top sub-domain, one third to the central one and one third to the bottom subdomain. At the beginning of coking there exist neither the top nor the bottom sub-domain. Using the numerical procedure described in Section 2 and the boundary conditions of the second kind (heat fluxes given) specified at the coal bed top and bottom, we calculate the temperature distribution at instance Dt. Then, the velocities of both the top and bottom evaporation fronts can be evaluated using Eqs. (28) and (29), respectively and the locations of the top and the bottom evaporation fronts can be evalu- Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 6 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 358 359 360 362 363 364 366 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx ated using Eq. (30). Now, using the first order spline we calculate (by linear interpolation) the temperatures at all the numerical elements (nodes) located in both the upper and bottom sub-domains. Then, the devolatilization sub-model (see below) is called to determine the composition of volatiles in each node (numerical cell) of the upper and bottom sub-domains. Subsequently, the solid bed properties (porosity, composition, true and apparent densities, specific heat, thermal conductivity and so on) are also updated in the top and bottom sub-domains. At the next instance, corresponding to 2Dt, the new location of the evaporation front is calculated and consequently the size of the sub-domain increases. Then, the heat transfer Eq. (1) is solved again in the top sub-domain using at its top the boundary condition of the second kind (specified heat flux) and at its bottom (at the evaporation front) the boiling temperature. Again, since the numerical grid at the instance 2Dt is different to the grid at instance Dt, the spline is needed for interpolation. In the central sub-domain the heat transfer equation (Eq. (1)) is solved using at the sub-domain boundaries a fixed temperature equals to Tboiling (boundary conditions of the first kind). The above described moving boundary technique is a bit complex and requires a tedious book-keeping while coding. Thus to summarize, Eq. (1) is solved in the three sub-domains in the following formulations: BOTTOM SUB-DOMAIN: for 0 < x < xBevap ðsÞ 8   d > qdbulk cdbulk @T@bulk ¼ @x@ keff;bulk @T@xbulk þ Sdev > s > d d < qð0; sÞ ¼ boundary condition at the bottom >   > > : T bulk xBevap ; s ¼ T boil   8 m > qmbulk cmbulk @T@bulk ¼ @x@m keff;bulk @T@xbulk þ Scond > s m > > <   T bulk xBevap ; s ¼ T boil > >   > > :T T bulk xevap ; s ¼ T boil ð33Þ TOP SUB-DOMAIN: for xTevap ðsÞ < x < H   8 d > qdbulk cdbulk @T@bulk ¼ @x@ keff;bulk @T@xbulk þ Sdev > s > d d <   T bulk xTevap ; s ¼ T boil > > > : qðH; sÞ ¼ boundary condition at the top Condensation occurs when the temperature of the coking-bed numerical cell is lower than Tboiling. The mass flow rate of moisture condensated is determined using the following formula: eVMm pboil sat pTsat
375 _ cond ¼ m 376 where e – coal/coke porosity (m3gas /m3bed ), V – cell volume (m3), Mm – 379 380 381 383 384 385 2.4. Devolatilization 399 Devolatilization occurs in the top and bottom sub-domains, as shown in Fig. 4. Thus, in every numerical cell of these domains the Distributed Activation Energy Model (DAEM), as proposed by Merrick [25], is used (see Figs. 5 and 6). As in Merrick [25], the composition of the volatile matter is defined in terms of the following nine species: CH4, C2H6, CO, CO2, tar, H2, H2O, NH3, H2S. For each species, the cumulative amount ‘‘mi” released at time ‘‘s” is given by: 400 Z s0 Z 1 j
exp E0  E RT bulk ðsÞ  i ½m mi ðE; s RT bulk Ds ð35Þ pboil sat pTsat molar mass of the moisture (kgm/kmol), – saturated vapour pressure at boiling temperature (N/m2 Þ, – saturated vapour pressure at cell temperature (N/m2). The saturated vapour pressure is given by Clausius–Clapeyron equation: pTsat   M m revap 1 ¼ patm exp R T boil 1 T bulk  ð36Þ where Tboil is boiling temperature (K) and pTsat is the saturation pressure (Pa). As shown in Fig. 4 there are two condensation fronts: 387 390 391 392 393 394 395 396 397 398 401 402 403 404 405 406 407 408 ފ
F 0i ðEÞ
dE
d s ð37Þ 410 where j – frequency factor (1/s), E – energy activation (J/kmol), R – universal gas constant (J/kmol K), Tbulk – deposit temperature (K), E0 – ‘‘starting” activation energy (kJ/kmol), F0 (E) – noncumulative Rosin–Rammler function, mi – final yields of coke and volatiles (kgi/kgdaf), mi (E) – cumulative volatile release at time t (kgi/kgdaf), s – total time (s). The non-cumulative Rosin–Rammler distribution has a form: 411 ð34Þ 2.2. Moisture condensation 378 389 F 0i ðEÞ ¼ 369 377 Formation of tars, its flow and condensation are important for the carbonization process. The whole process of tar formation and condensation is much more complicated than the above described transport of water. The complexity stems from variable composition of tars. In this work the tars are considered by the use of additional energetic effects (heat sources and sinks) obtained using the inverse procedure described in Part IV (see [67]). Indeed, this is a simplification but perhaps it is fair to say that in order to model tar formation and condensation a separate, very complex, research program is required. CENTRAL SUB-DOMAIN: for xBevap ðsÞ < x < xTevap ðsÞ scripts d; m indicate the dry and wet coal, respectively. 372 388 0 368 373 2.3. Tar formation and condensation mi ðE; sÞ ¼ where xTevap ; xBevap are positions of the evaporation fronts, the sub- 371 386 ð32Þ 367 370 _ Tcond (moving downwards) and the second at the one on the top m B _ cond (moving upwards). bottom m b E e E0  e b 1 exp (  E E0 e b ) ð38Þ where E – activation energy (J/kmol), e ¼ e1 e2 – width parameter, b – shape parameter, E0 – ‘‘starting” activation energy (J/kmol). In general, the parameters E0, e, b vary with gas components. The values used in the present work are shown in Table 1. The double integral (37) has to be calculated in every numerical cell of the upper and bottom sub-domains at each time step Ds, see Figs. 5 and 6. The integral can be simplified to the form: Dmi ðEÞ ¼ Z Ds Z 0 1 j
exp E0  h E i m RT bulk ðsÞ mi i
F 0i ðEÞ 412 413 414 415 416 417 418 420 421 422 423 424 425 426 427 428
dE
ds ð39Þ 430 where E – activation energy (J/kmol), R – universal gas constant (J/kmol K), Tbulk – the bulk cell temperature (K), E0 – ‘‘starting” activation energy (J/kmol), F0 (E) – non-cumulative Rosin–Rammler  i – final yields of coke and volatiles (kgi/kgdaf), function, m 431 ðsÞ mi – cumulative volatile release at time s (kgi/kgdaf). Performing the integration over residence time of the solid matter in the numerical cell, the double integral becomes a single one: Dmi ðEÞ ¼ Z 1 E0 j
exp h E i m RT bulk  ðsÞ mi i
F 0i ðEÞ
dE
Ds ð40Þ In every time step and in each cell for each component ‘‘i” the cumuðsÞ lative volatile release mi is calculated and stored in memory: Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 432 433 434 435 436 437 438 440 441 442 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 7 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 5. Usage of Merrick’s devolatilization model [25] (top of the computational domain). 443 ðsþDsÞ ðsÞ 445 mi 446 where ‘‘i” represents the i-th volatile component. Composition (mass fractions) of the gas during devolatilization at certain time s can be evaluated using: 447 448 ¼ mi þ Dmi ðEÞ ð41Þ where Dmi is related to mass of the dry ash free coal and expressed in kgraw/kgdaf. The subscript dev indicates raw pyrolysis gas. The amount of the raw-gas produced in the coking-bed is obtained upon integration over the top and bottom sub-domains, so that 452 453 Dmi ðEÞ g dev ¼ P9 i i¼1 Dmi ðEÞ ð42Þ The final yield (mi) of each species and coke (see Table 2) are calculated using the following system of equations [25]: 472 473 474 475 449 451 471 mdev n X ¼ mk;dev
qk;true
ð1 ek Þ
ð1 ak wk Þ
V k ð45Þ k¼1 where mdev is expressed in (kg), qtrue is the true density of coal/coke (kg/m3), e is the fixed bed porosity (m3/m3), a, w – mass fractions of 3 32  3 2 mchar c 0:98 0:75 0:8 0:4286 0:2727 0:85 0 0 0 0 7 6 7 7 6 0:002 0:25 0:2 6  CH4 7 6 h 0 0 0:082 1 0:1111 0:1765 0:0588 76 m 7 6 7 76  7 6 6 7 6 0:002 0 0 0:5714 0:7273 0:049 0 0:8889 0 0 76 mC2 H6 7 6 o 7 76 7 6 6 7 6 76 m 7 6 0:01  n 0 0 0 0 0:009 0 0 0:8235 0 7 76 CO 7 6 6 7 7 6 76 6  7 6 7 7 6 0:006 6 mCO2 s 0 0 0 0 0:01 0 0 0 0:9412 7 7¼6 76 6 6 7 7 6 1 6  tar 7 6 1 V 7 0 0 0 0 0 0 0 0 0 76 m 7 6 7 7 6 76 6  H2 7 6 1:31h 7 76 m 6 0 1 0 0 0 0 0 0 0 0 7 7 6 76 6 7 6 0:22h 7 6 0 6  0 1 0 0 0 0 0 0 0 7 H2 O 7 7 6 6 76 m 7 7 6 6 76 4 0  NH3 5 4 0:32o 5 0 0 1 0 0 0 0 0 0 54 m  H2 S m 0:15o 0 0 0 0 1 0 0 0 0 0 477 478 479 2 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 470 where c, h, o, n, s, V – mass fractions of carbon, hydrogen, oxygen, nitrogen, sulphur and volatile matter in the daf coal, V = p 0.36p2 (correction introduced by Merrick [25]), where p is the ASTM volatile matter content. The first five equations are the balances of carbon, hydrogen, oxygen, nitrogen and sulphur. The sixth equation represents the overall mass balance. The remaining four equations are associated with additional inter-correlations described in [25]. In order to solve the entire system of equations one has to know the proximate and ultimate analysis of the coal. Moreover, the composition of the final coke and tar is required and these are given in Table 2. The amount of (raw) gas released during devolatilization process from the single computational cell (k), at a given time, is calculated as follows: mk;dev 9 X ¼ Dmi ðEÞ ð44Þ ð43Þ ash and moisture respectively (kg/kg), V is the cell volume (m3). The summation extends over all numerical cells. The raw-gas flow rate, at a given time, is calculated as follows: _ dev ¼ m mdev Ds ð46Þ One should bear in mind that two devolatilization zones exist: _ Tdev and the bottom one generating m _ Bdev . the top one producing m Moisture evaporation occurs at the two internal moving boundaries (evaporation fronts, see Fig. 4) whereas moisture condensation takes place in the central sub-domain. To determine the raw-gas flow rate leaving the coking-bed (input to the sub-model describing the space above the coking-bed, see Part III [32]) one adds up the streams due to evaporation, condensation and devolatilization. _ raw ¼ m _ Tevap m þ _ Bevap m _ TC m cond _ BC m cond þ _ Tdev m þ _ Bdev m ð47Þ i¼1 Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 480 481 482 483 485 486 487 488 489 490 491 492 493 494 495 497 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 8 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 6. Usage of Merrick’s devolatilization model [25] (bottom of the computational domain). Table 1 The parameters adopted in the Rosin–Rammler distribution [25]. Eo (MJ/kmol) b e1 (MJ/kmol) e2 (MJ/kmol) CH4 C2H6 CO CO2 Tar H2 H2O NH3 H2S 183 2 110 0 183 4 61 0 183 4 93 0 183 4 78 0 183 8 23.6 17.6 183 4 165 0 183 8 23.6 17.6 183 4 106 0 183 4 114 0 Table 2 Ultimate analysis of coke and tar [25]. Tar Coke Carbon Hydrogen Oxygen Nitrogen Sulphur 0.85 0.98 0.082 0.002 0.049 0.002 0.009 0.01 0.01 0.006 504 where T denotes the top (crown side), B the bottom (sole-flue side), whilst TC indicates the moisture evaporated in the top sub-domain and condensated in the central sub-domain; similarly BC indicates the moisture that evaporated in the bottom sub-domain and condensated in the central sub-domain. In addition, the composition (mass fractions) of the raw-gas is calculated as 507 g raw ¼ g dev i i 508 and for H2O: 498 499 500 501 502 503 505 509 _ dev m _ raw m dev;T g raw H2 O ¼ g H2 O 513 514 515 517 518 519 520 521 522 523 524 The coal/coke composition has an influence on physic-chemical properties of the charge. Properties such as thermal conductivity, specific heat, true density (see Section 3) and porosity are necessary to correctly determine the temperature profile. Mineral matter is a part of the coal/coke charge and hence it has been included in the model as ash component. It was assumed that the ash is inert and therefore is not released into gas phase. 525 2.6. Coking-bed composition 532 During the unsteady calculations of the carbonization process it is necessary to store the information concerning the composition of the coking-bed. In each numerical cell of the bed the ash, water and combustible part (char plus volatiles remaining in the char are regarded as coke) content have to be known. The mass fractions of moisture xm, ash xa and combustibles (coke) xcp in the coking-bed are calculated using the following formulae (we annotate mass fractions of the coking-bed by x to distinguish from mass fraction in gas phase annotated by g, see the previous section): 533 xðmsþDsÞ ¼ ð48Þ q q
1 e ðsÞ 526 527 528 529 530 531 534 535 536 537 538 539 540 541 542
ðsÞ V eðsþDsÞ ÞV ðsþDsÞ
_ cond Ds _ evap m m sþDsÞ qðtrue
ð1 eðsþDsÞ ÞV ðsþDsÞ ð51Þ 544 545 _ Tevap m _ Bevap _ Tdev _ Bdev m m m þ g Hdev;B þ þ O _ 2 _ _ _ raw mraw mraw mraw m _ TC m cond _ raw m xðasþDsÞ _ BC m cond _ raw m ¼ where gi is the mass fraction of the i-th raw-gas component. The average raw-gas temperature can be estimated from the following formula: Pn k¼1 mk;dev c k;bulk  T k;bulk mdev cbulk To
ð50Þ where n – number of grid elements, ck;bulk – mean specific heat of the coal located in the k-th numerical cell (J/kg K), cbulk – mean specific heat of the entire coal bed (J/kg K), T k;bulk is the coal temperature in the k-th cell (K), T o = 298.15 K is the reference temperature. The composition, temperature and amount of raw-gas are then used in the sub-model describing the upper-oven (see Part III [32]). ðsÞ xa
qtrue
1 sþDsÞ qðtrue
ð1 eðsÞ V ðsÞ
eðsþDsÞ ÞV ðsþDsÞ ð52Þ 547 548 ð49Þ T raw ¼ T o þ ðsÞ ðsÞ xm
true ðsþDsÞ
ð1 true ðsÞ 511 512 2.5. Mineral matter ðsÞ ðsþDsÞ xcp ¼ ðsÞ xcp
qtrue
1 sþDsÞ qðtrue
ð1 eðsÞ V ðsÞ
eðsþDsÞ ÞV ðsþDsÞ _ dev Ds m ðsþDsÞ qtrue
ð1 eðsþDsÞ ÞV ðsþDsÞ ð53Þ where x – mass fraction of the bed component (kg/kg), m – moisture, a – ash, cp – combustible part (coke), ðs þ DsÞ – next time step, ðsÞ – current time step. The model assumes that during non-stationary calculations the volume of numerical cells is fixed; it means the bed does not shrink. The porosity of each numerical cell is calculated for every time step: Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 550 551 552 553 554 555 556 557 558 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 560 sÞ qðtrue
1 ðsþDsÞ ðsþDsÞ qtrue V
eðsÞ VðsÞ eðsþDsÞ ¼ 1 s _ raw D m ðsþDsÞ ðsþDsÞ true V q ð54Þ 561 2.7. Coke composition 562 Composition of the coke (char plus volatiles in the char) can be determined using the existing system of equations (see Eq. (43)). The system (43) can be written as follows [25]: 563 564 565 567 568 569 570 571 572 9 X  j ¼ bi ; Aij m 577 578 580 581 582 584 9 X mj ðsÞ ð56Þ j¼1 where index j indicates the j-th volatile matter component. The composition of the char residues at time ‘‘s” is calculated by element balances (the first five equations of system (43) [25]): 9 X yi m0 ðsÞ þ Aij mj ðsÞ ¼ bi ; i ¼ 1; . . . ; 5 ð57Þ j¼1 that is: yi ¼ bi P9 j¼1 Aij mj ð m0 ðsÞ sÞ ; i ¼ 1; . . . ; 5 ð58Þ 587 Thus, mass fractions of c, h, o, n, s in the coke and mass fraction of volatiles contained in the coke are calculated in each numerical cell, at each instant of the carbonization process. 588 3. Thermal properties of the coking-bed 589 In order to describe the heat transfer process in the coking-bed, one has to take account of substantial variations in thermal properties with coking time. The specific heat (cbulk Þ and the thermal conductivity (keff;bulk Þ strongly depend on temperature of the 585 586 590 591 592 593 3.1. Thermal conductivity 603 Thermal conductivity is one of the most important properties in modelling of fixed-beds. Fig. 7 shows thermal conductivity values which originate from different literature sources. For hard coals (see Fig. 7 curves D–E) the conductivity at 25 °C ranges from 0.1 to 0.5 W/mK (depending on the coal grade) and at 1000 °C thermal conductivity of coke produced increases only slightly to 0.2–0.6 W/mK. Coal beds are good insulators and reluctantly convey heat. Resistances to heat transfer occur at contact surfaces between solid and gaseous phases as well as in pores and fissures (cracks) between grains that trap air (thermal conductivity around 0.025 W/mK). When the fixed-bed temperature exceeds 600 °C, the radiation between particles prevails, which increases the heat transfer rate several times. Fig. 7 shows the effect of radiation on the so-called effective thermal conductivity (thermal conductivity + effect of radiation) (see Fig. 7 curves A–C and F–G). The effective thermal conductivity increases significantly with temperature; at 900 °C its values are in the 0.8–2.0 W/mK range. Fig. 7 also shows large discrepancies in the thermal conductivity values quoted in different literature sources. Significant discrepancies appear above 600 °C, in a temperature region where thermal radiation is dominant. Thus, establishing a relationship for calculating the effective thermal conductivity is a sophisticated problem since the bed structure varies with time. It is difficult to estimate shapes and dimensions of pores where the radiation process takes place. Many publications simplify the sophisticated structure of the fixed-bed and restrict the analysis to idealized geometrical forms made up of spherical particles. The simplification may lead to conductivity values departing substantially from reality [28]. It has been established that the effective thermal conductivity depends on such parameters as [35–37]: (a) particle diameter and emissivity, (b) dimensions, size and type of pores, (c) porosity of the fixed-bed, (d) temperature of the bed. Hütter and Kömle [36] have established that the effect of radiation is substantial when the temperature of coal/coke-bed exceeds 600 °C and particles are larger than 1 mm. Computations of Schotte [37] examine the effect of pore sizes and porosity of the fixed-bed on the value of the radiative part of thermal conductivity. The work demonstrates that larger contact surfaces between the solid and gaseous phases lead to more intense transfer of heat via radiation. In fact it is the case when porosity of the bed is high and the pores are small. The research report of Atkinson and Merrick [28] describes in a very accurate and detailed manner how heat is transferred in coal/ coke-beds. After determination of the effective thermal conductivity, that takes account of convective heat exchange in pores as well as heat transfer by radiation, attention is paid to amendments in coal/coke structure and composition during the combustion process. Three basic forms of solid-beds have been distinguished. For the first, original form, one assumes that particles are not porous which means that the fixed-bed contains only external pores. The second and third forms (plastic and sintered forms) not only comprise external pores but also internal ones that are created during 604 594 595 596 597 598 599 600 601 602 ð55Þ where A – matrix of constants, b – vector of constants, m – vector of final yields of coke and volatiles. The cumulative masses of the volatile matter components are used to determine the mass of char moðsÞ remaining at time s [25]: 574 576 material as well as on the bed composition and structure. Similarly to thermal conductivity and specific heat, also bulk density (qbulk ) is subjected to variations. The bulk density of the charge decreases with increasing temperature since the total porosity (e) of the bed rises with temperature. In our model calculations true density ðqtrue Þ of the solid matter of the coking-bed increases with the coking time because species of relatively low density, such as moisture and volatiles are released to gas phase (vaporization, devolatilization). In reality, the true density changes are also associated with rearrangement of the coal organic matter. j¼0 m0 ðsÞ ¼ 1 575 i ¼ 1; . . . ; 10 9 Fig. 7. Thermal conductivity of coal/coke beds. C, D, E – thermal conductivity (compiled by Tomeczek [33]), B, F – effective thermal conductivity (compiled by Tomeczek [33]), A – effective thermal conductivity [28], G – effective thermal conductivity used in this work [34]. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 10 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx estimated in several literature sources [26,34,38–44], for a number of coals. The relationships differ from each other; not only shapes of the curves are different but also the maximum occurs at different temperatures. Observed discrepancies are mainly due to the fact that different coals are considered. In some works, as for example in publication of Eisermann et al. [38] (curve I in Fig. 8), the exothermic/endothermic effects of carbonization have been included. Curve ‘‘I” in Fig. 8 shows that the endothermic effects in the temperature range 150–800 °C have been accounted for by an increased value of the effective specific heat. Curve ‘‘A” in Fig. 8, which originates from the report of Atkinson and Merrick [26], shows the effect of temperature on the equivalent specific heat of the solid matter (coal/coke) forming the bed. The measurements [26] have demonstrated that the specific heat of the coal/coke considered reaches a maximum of 2200 J/kg K at temperatures around 500 °C. A further temperature increase leads to a drop which is caused by releasing substances with the highest values of specific heat, such as volatiles. Atkinson and Merrick have observed [26] that the endothermic effects may actually not occur at temperatures below 700 °C. Tomeczek [33] proposes an experimental relationship that determines the effect of the temperature on the specific heat for moisture, volatiles, char and ash. The suggested method reproduces actual variations of thermal capacity during combustion of solid fuel and confirms observations reported in [26]. In this study an approach proposed by Kim in [34] is used (see Fig. 8, curve J). The function used does not differ from the other dependences to a large extent. The values of specific heat vary from 1.1 kJ/kg K at the beginning of the process to 1.8 kJ/kg K at the end of the process. Its maximum value (2.2 kJ/kg K) is reached at temperature of about 750 K. 698 687 physical transformation of coal as well as during coal devolatilization. After final extraction of tar and volatiles the bed achieves its final form. The remaining char cracks, due to thermal effects, which leads to formation of pores of very specific shapes. Atkinson and Merrick [28] have distinguished two major zones to enable mathematical description of heat transfer by conduction and radiation. In the first one, named as particulate charge, the external pores exist only. In the second zone, named as coke charge, both types of porosity exist; the original fixed bed porosity and the extra porosity due to the cracks formed. Heat transfer proceeds in these two zones in different ways. In the first zone, thermal conduction in coal and moisture, conduction in pores and radiation between fuel particles take place. Within the second (coke) zone the effective thermal conductivity has to be determined using a different procedure due to the fact that two types of porosity (internal and external ones) exist. The pores have different shapes as compared to the space between particles. Alteration of both shapes and dimensions of pores affect the characteristic dimension that is used to determine the radiative part of thermal conductivity. The coke zone features heat conduction in coke or gas, radiation across internal pores as well as along cracks. The approach proposed in [28] describes the heat transfer within the fixed-bed in a strict and accurate manner, as the authors adopted the model that takes account of such important phenomena as alteration of solid fuel structure during the combustion process. The dependence proposed in [34] is used in this study (see Fig. 7, curve G) to calculate the effective thermal conductivity. At temperatures lower than around 800 K the conductivity, calculated using dependence G, is lower if compared to the other relationships. It is larger at temperatures larger than around 800 K. Thus, in our work, the influence of radiative heat transfer on the effective thermal conductivity is larger than in the cited literature sources. 688 3.2. Specific heat 3.3. True density 729 689 Another important property is the specific heat of the bed. Composition of the coking-bed is subject to substantial variations due to moisture evaporation and condensation, pyrolysis and carbonization. Each component (moisture, volatiles, char and ash) possesses a specific heat that is different. It leads to variations of specific heat with fuel composition. Vast majority of published works uses the equivalent specific heat where the bed is considered as the mixture of moisture, volatiles, char and ash [26,33]. Fig. 8 shows the effect of temperature on equivalent specific heat Similarly to thermal conductivity and specific heat, also true density of solid matter is subject to variations during carbonization. The increase of the bed temperature results in an increase of bed density as chemical compounds with relatively low density, such as water and volatiles, are released to gas phase (vaporization, devolatilization). The solid phase contains then more char and ash. The relationship between the coal/coke true density and temperature is shown in Fig. 9. Fig. 9 shows clearly that the true density of four different coals does not differ significantly from each other although differences can be seen at temperatures below 800 °C 730 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 690 691 692 693 694 695 696 697 Fig. 8. Effective specific heat of fixed coal/coke beds (A – [26], B, C – [39,40], D – [41], E – [42], F – [43], G – [44], H, I – [38], J – relationship used in this work [34]. Fig. 9. Properties of the coal/coke bed – true density. B, C, D – [27], A – relationship used in this work [34]. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 731 732 733 734 735 736 737 738 739 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 11 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 10. Coking-bed porosity. A – conventional slot-type oven [27], B – HR/NR horizontal oven [34]. 745 (compare curves C, D and B). A high volatile coal has a lower density (see curve D) if compared to that of a low volatile coal (see curve C) [27,34]. The relationship proposed by Kim in [34] is used in this study (see Fig. 9, curve A) and the density values are similar to those reported in the literature (see Fig. 9). 746 3.4. Porosity 747 Fig. 10 shows changes to the bed porosity during carbonization. The overall volume of the pores with respect to the total volume of the bed increases with temperature from around 35% to around 70% (for conventional coke oven – see Fig. 10, curve A) and from around 12% to around 65% (for HR/NR coke oven – see Fig. 10, curve B). Total porosity of the coal/coke bed is lower in the case of HR/NR coke making process since the coal bed is highly pressed/ compacted. In the coking-bed sub-model formulated in this work, the porosity of the bed is calculated in each numerical cell and at each instance, see Eq. (54). However, an initial (at the beginning of the process) porosity is required which has been estimated to be 740 741 742 743 744 748 749 750 751 752 753 754 755 756 757 758 0.15. Such a low porosity (typical values are in the 0.3–0.4 range) is applicable here since the coal-bed (charge) has been pressed before coking. 759 3.5. Bulk density 762 Knowing both total porosity and true density, the bulk density can be calculated. Bulk density, in contrast to the true density, decreases during coke making process since the total porosity of the coal/coke bed rises. Charges of conventional coke ovens are characterized by bulk density in the range from around 900 kg/ m3 to around 600 kg/m3, as shown in Fig. 11 – curve A. Bulk density of the charge prepared for HR/NR coke ovens is higher and varies from 1150 kg/m3 to 750 kg/m3 (see Fig. 11, curve B). It is worth noting that in the400–550 °C temperature range (typical for plastic stage) the bulk density of both coals may temporarily increase [27,34]. In our modelling work the bulk density of the cokingbed is calculated using the porosity and true density. 763 3.6. Thermal effects of carbonization 775 Thermal effects accompanying carbonization have an impact on heating of the charge and they are coal dependent as well as heating rate dependent. Moreover, air availability may play a role. The analysis provided in Merrick [26] shows that thermal effects may be of two types. The endothermic effects of around 200 kJ per kg of charge occur at temperatures up to 900 K and can be neglected only when the oxygen content is sufficiently high. The exothermic effects of around 200–400 kJ per kg of charge occur at temperatures above 900 K. It has also been reported [26] that a rapid heating of the charge brings about (300–1400) kJ/kg energy released in exothermic reactions. Rycombel and Długosz [45] have investigated the thermal effects of devolatilization process of six powdered and compacted coals. The tests have been conducted for three different heating rates 3, 6 and 10 K/min under protective atmosphere (nitrogen, argon). The results show that the endothermic and exothermic effects are different for each coal type. Furthermore, thermal effects may occur in various temperature ranges. For most coals exothermic effects exist in the 160–400 °C temperature range and are at the level of 24–39 kJ/kg. However, in the 160–750 °C temperature range, the endothermic effects are dominant and they are at the level of 160–293 kJ/kg. It should be emphasized that these results are only appropriate for certain coals. In some cases, endothermic and exothermic effects can appear at different temperature ranges. 776 Table 3 Coking pressure for different bulk densities (conventional oven) [65]. Bulk density, qbulk (kg/m3) 720–750 820–860 920–950 1000 Coking pressure pcoking (kPa) Experiment According to Eq. (59) 15–18 20–25 34–37 15.8 21.5 34.5 47.7 Table 4 Plastic layer thickness and shrinkage for different coals [65]. Fig. 11. Bulk density of coking-bed; A – conventional slot-type ovens [27], B – NR/ NR horizontal ovens [34]. Yield of volatiles (%) Plastic layer thickness (mm) Gap at wall (mm) 22 31 27–29 13 12 17–19 3–4 6–9 20 Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 760 761 764 765 766 767 768 769 770 771 772 773 774 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 12 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 12. Coking pressure curves for several American coals (movable-wall test oven) [61]. 801 802 803 804 805 In the coking model presented in this paper, it has been assumed that during carbonization process, in 160–550 °C temperature range, the endothermic effects exist at the level of 150 kJ/kg. The exothermic effects appear at temperatures larger than 550 °C and are equal to 300 kJ/kg. Fig. 14. Peaks of internal gas pressures measured in selected locations of the charge (movable-wall test oven) [58]. 4. Coking pressure in HR/NR coke oven 806 At around 350 °C temperature, coal is subjected to depolymerization and transformation into the unsteady phase, becoming a sticky liquid. That phase, referred to as the softening phase, may 807 Fig. 13. Coking pressure and internal gas pressure (movable-wall test oven) [67]. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 808 809 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 15. Coking sub-model – computed temperature profile along charge height. 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828 829 830 831 832 833 834 835 836 837 838 839 occur or not, depending on the fuel grade, heating rates, dimensions of particles and composition of the gaseous atmosphere where the fuel is placed. The plastic (soft) phase is associated with the phenomenon of coal particle swelling during heating up at high rates. The plastic layer formed in the charge affects behaviour of water vapour and volatile matter within the bed. According to Foxwell [46] the plastic layer creates so-called plastic envelope which is partially impermeable to gas flow so that the coking pressure increases. Degree of impermeability depends on strength (viscosity) and thickness of the plastic layer. In other publications [47– 49] it is argued that the release of moisture and volatiles to the space above the coking-bed is restricted not by plastic layer but by coke and semi-coke. Very popular is the view that the pressure is generated directly in the plastic layer [50–54] while others assume that the pressure occurs inside swelling coal particles when pyrolysis occurs [55–57]. The impermeable plastic layer blocks evaporated water and volatiles which then condensate in a cold region in the centre of coal bed (inside the envelope) increasing the coking pressure. Several investigators [58,59] noticed that the internal gas pressure inside the envelope rises at about 100 °C and depends on (a) coalification of the coal, (b) coal coking properties, (c) particle size, (d) bulk density and (e) coking rate. Coking pressure has been measured to be in the 1–17 kPa range for full-scale industrial ovens [60–62] and up to 200 kPa for smallscale, laboratory, test installations [63]. The pressure inside the closed pores of coal particles may be higher than the pressure in the bed centre [64,65]. Gagarin [65] provides dependence of the mean coking pressure on the bulk density of the coal charge (see Table 3) in the form 13 shown as wall pressure, increases with coking time whilst the internal gas pressure rises only for a moment, when plastic layer pass through the measuring probe. The internal gas pressure is higher than the pressure monitored at the oven walls. A similar pressure peak may be found during evaporation of moisture [67]. Fig. 14 shows two pressure curves; one for the oven wall pressure and another for the internal gas pressure. The latter has been measured at five locations (marked as P1–P5) inside the charge. Fig. 14 indicates that internal gas pressure rises significantly in the moving plastic layer. When the plastic layer passes the pressure probe (sensor) there is a sudden rise in the gas pressure. As soon as the plastic layer moves beyond the probe, the gas pressure drops quickly and a rise in the gas pressure is not recorded anymore. The maximum of the internal gas pressure (35 kPa) is achieved at the centre of the bed. As in South and Russel [75], the bulk bed density of 750–850 kg/m3 [58] is applicable for curves shown in Fig. 14. Measurements in HR/NR coke ovens have also shown that the internal pressure raises during carbonization and it affects the moisture boiling temperature. In the coking-bed sub-model presented in this paper the boiling temperature increases linearly from 373 K at the beginning of the coke making process to the 383 K at the end, which corresponds to water saturation pressures of 101 kPa and 142.5 kPa, respectively (see Eq. (36)). 856 5. Numerical simulations of the carbonization process 880 The above described coking-bed sub-model is used together with the hydraulic network sub-model [1] and the other submodels [32] handling the upper-oven, the down-comers and the sole-flues. This facilitates time-dependent simulations of the oven performance and the coking-bed sub-model is being run at continuously varying boundary conditions. In other words, heat fluxes at the top and at the bottom of the coking-bed vary with time. Results of such simulations are presented in Part IV [71]. In this section we demonstrate the capabilities of the coking-bed sub-model only. To this end, the coking-bed sub-model is being run assuming constant temperatures at the bottom and at the top of a 1 m high bed; both temperatures are taken to be 1200 K. The initial coal-bed temperature equals 298 K, the initial true coal density is 1300 kg/m3, whilst a value of 1.1 kJ/kg K is used for the initial coal specific heat. The initial bed thermal conductivity is taken to be 0.15 W/mK whilst the initial porosity is 0.15. As described above, all these 881 840 0:5942qbulk þ 0:000412q2bulk 842 pcoking ¼ 229:921 843 where the bulk density is in kg/m3 while the pressure in kPa. Table 4 shows that plastic layer thickness raises with the volatiles yield. The measurements indicate shrinkage of coke cake while heating which in Table 4 is shown as thickness of the gap formed between the oven walls and the bed. The shrinkage increases with volatiles yield. Soth and Russel [61], see Fig. 12, have presented coking pressure measurements in the centre of the coal/coke charge for several American coals crushed below 3 mm and charged with bulk density of 785–850 kg/m3. When the plastic layer approaches the centre of the charge the coking pressure approaches a maximum at 35 kPa. Similar results have been obtained at the BCURA [66]. Fig. 13 demonstrates a wall pressure and a gas pressure measured in the charge centre for three coals. The overall coking pressure, 844 845 846 847 848 849 850 851 852 853 854 855 ð59Þ Fig. 16. Calculated composition of the raw gas produced during carbonization. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 857 858 859 860 861 862 863 864 865 866 867 868 869 870 871 872 873 874 875 876 877 878 879 882 883 884 885 886 887 888 889 890 891 892 893 894 895 896 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 14 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx Fig. 17. Mass flow rate of the raw gas released to the space above the coking-bed (half of the oven). properties are recalculated, in each numerical cell, at any instance of coking time. The thermal effects associated with coking are included, as described in Section 3.6. Fig. 15 shows calculated temperature profiles along the bed height; the curves correspond to different coking times. As can be seen, at the beginning, the coal/coke bed is heated up to the water boiling temperature which varies due to an increase in the internal pressure of the charge. The evaporation fronts move from the bottom (sole-flue side) and the top (crown side) towards the charge centre with 0.01 mm/s and 0.003 mm/s velocities at the beginning and at the end of the evaporation process, respectively. A certain amount of moisture condensates in the (cold) part of the bed. Fig. 15 shows that after moisture evaporation a rapid heating occurs at the top and at the bottom of the bed. The evaporation fronts meet in the centre of the charge after about 56 h and this instance corresponds to complete evaporation of moisture. As soon as the whole moisture is given off the bed temperature raises rapidly. The middle of the bed reaches the temperature of 1200 K after about 80 h. Fig. 16 shows the calculated composition of the raw-gas. At the beginning of the carbonization a large amount of moisture is Fig. 18. Coking-bed composition (coke, ash, moisture) at three instances. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 897 898 899 900 901 902 903 904 905 906 907 908 909 910 911 912 913 914 915 916 917 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 15 Fig. 19. Coke elemental composition during carbonization. 918 919 920 921 922 923 924 925 926 927 928 929 930 931 932 933 934 935 936 937 released. After about 10 h the gas composition stabilizes until the evaporation fronts meet in the centre of the bed (56 h). At this instance the moisture has been completely removed and a rapid increase of the bed heating rate is observed. The amount of rawgas produced (calculated for half of the oven) also rapidly increases at 57 h instance, as shown in Fig. 17. In the last ten hours of the carbonization, mainly hydrogen is released albeit in very small amounts. Fig. 18 shows the calculated composition of the coking-bed at three selected instances in coking time. More precisely, the mass fraction of coke (char plus the volatiles remaining in the char), ash and moisture are plotted. It can be seen that a certain portion of evaporated water condensates in the cold areas which increases the humidity of the solid phase and decreases the evaporation rate. Fig. 18 shows that after evaporation of moisture the volatile matter is quickly released. After carbonization is completed the bed contains 90% coke and 10% ash, by weight. Fig. 19 shows the elemental composition of the coke during carbonization. It may be noted that almost all hydrogen and oxygen are released into the gas phase. At the end of coking, the coke con- tains 98% of carbon, 2% of hydrogen and 2% of oxygen, by weight as specified in Table 2. 938 6. Final remarks 940 This paper is a part of the series of articles ‘‘One-dimensional model of heat-recovery, non-recovery coke ovens” [1,32,67]. In this paper we describe the development of a one-dimensional timedependent sub-model of the coking-bed. The sub-model predicts the temperature distribution, composition (mass fractions of coke, ash and moisture) as well as physical structure (porosity, density) of the coking-bed as a function of time and position in the bed. In order to handle the discontinuities occurring due to moisture evaporation and condensation we use a moving boundary method to solve the heat balance equation and this is indeed a novel approach. The coal pyrolysis has been modelled using Merrick et al. [25] approach. The coking-bed sub-model provides qualitatively correct predictions for 80 h coking process of a compact charge representative of those typical for NR/NR ovens. Detailed 941 Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 939 942 943 944 945 946 947 948 949 950 951 952 953 954 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 16 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 956 validation and assessment of the overall HR/NR model, and the coking sub-model in particular, are provided in Part IV [67]. 957 Acknowledgements 958 967 Since 2006 TU Clausthal and the Division of Coke Plant Technologies of ThyssenKrupp Industrial Solutions AG (TK IS) have been developing heat recovery coke ovens. After completion (in 2010) of the first phase concerning optimization of the primary air inlets within the upper-oven, the work has focused on the mathematical model development with the aim of optimizing the sole-flues design and the secondary air inlets. This paper is a part of this development. We would like to express our sincere thanks to ThyssenKrupp Companhia Siderúrgica do Atlântico, Brazil for their generous support. 968 References 955 959 960 961 962 963 964 965 966 969 970 971 972 973 974 975 976 977 978 979 980 981 982 983 984 985 986 987 988 989 990 991 992 993 994 995 996 997 998 999 1000 1001 1002 1003 1004 1005 1006 1007 1008 1009 1010 1011 1012 1013 1014 1015 1016 1017 1018 1019 1020 1021 1022 1023 1024 1025 1026 1027 1028 [1] Buczynski R, Weber R, Kim R, Schwöppe P. One-dimensional model of heatrecovery, non-recovery coke ovens. Part I: General description and hydraulic network sub-model. Fuel; 2016 [submitted for publicaion]. [2] Koizumi M. The combustion of solid fuels in fixed beds. The sixth symposium (international) on combustion, vol. 6, no. 1. The Combustion Institute; 1957. p. 577–83. [3] Cho Y, Joseph B. Heterogeneous model for moving-bed coal gasification reactors. Ind Eng Chem Process Des Dev 1981;20:314–8. [4] Radulovic P, Usman Ghani M, Smoot L. An improved model for fixed bed coal combustion and gasification. Fuel 1995;74(4):582–94. [5] Cooper J, Hallett W. A numerical model for packed-bed combustion of char particles. Chem Eng Sci 2000;55:4451–60. [6] Sudiro M, Pellizzaro M, Bezzo F, Bertucco A. Simulated moving bed technology applied to coal gasification. Chem Eng Res Des 2010;88:465–75. [7] He C, Feng X, Khim HC, Li A, Liu Y. Industrial-scale fixed-bed coal gasification: modeling, simulation and thermodynamic analysis. Chin J Chem Eng 2014;22 (5):522–30. [8] Buczyński R, Weber R, Szlek A, Nosek R. Time-dependent combustion of coal in a fixed-bed: measurements and mathematical modelling. Energy Fuels 2012;26:4767–74. [9] Ross AB, Jones JM, Chaiklangmuang S, Pourkashanian M, Williams A, Kubica K, et al. Measurement and prediction of the emission of pollutants from the combustion of coal and biomass in a fixed bed furnace. Fuel 2002;81 (5):571–82. [10] Oman J, Tacer M, Tuma M. Overfeed fixed-bed combustion of wood. Bioresour Technol 1999;67:139–47. [11] Thunman H, Leckner B. Modeling of the combustion front in a countercurrent fuel converter. Proc Combust Inst 2002;29:511–8. [12] Thunman H, Leckner B. Co-current and counter-current fixed bed combustion of biofuel a comparison. Fuel 2003;82:275–83. [13] Bruch C, Peters B, Nussbaumer T. Modelling wood combustion under fixed-bed conditions. Fuel 2003;82(6):729–38. [14] Zhou H, Jensen A, Glarborg P, Jensen P, Kavaliauskas A. Numerical modeling of straw combustion in a fixed bed. Fuel 2005;84:389–403. [15] Kolb T, Bleckwehl S, Gehrmann H, Seifert H. Characterisation of combustion behaviour of refuse derived fuel. J Energy Inst 2008;81(1):1–6. [16] Scharler R, Benesch C, Neudeck A, Obernberger I. CFD based design and optimization of wood log fired stoves. In: De Santi GF, editor. 17th European biomass conference & exhibition, June 2009. Florence, Italy; Hamburg, Germany: ETA-Renewable Energies; 2009. p. 1361–7. [17] Kausley S, Pandit A. Modeling of solid fuel stoves. Fuel 2010;89(3):782–91. [18] Collazo J, Porteiro J, Míguez JL, Granada E, Gómez MA. Numerical simulation of a small-scale biomass boiler. Energy Convers Manage 2012;64:87–96. [19] Buczyński R, Weber R, Szle˛k A. Innovative design solutions for small-scale domestic boilers; combustion improvements using a CFD-based mathematical model. J Energy Inst 2015;88:53–63. [20] Ryfa A, Buczynski R, Chabinski M, Szlek A, Bialecki R. Decoupled numerical simulation of a solid fuel fired retort boiler. Appl Therm Eng 2014;73 (1):794–804. [21] Athanasios N, Athanasios N, Nikolaos N, Nikolaos M, Panagiotis G, Kakaras E. Optimization of a log wood boiler through CFD simulation methods. Fuel Process Technol 2015;137:75–92. [22] Kær S. Numerical modeling of a straw-fired grate boiler. Fuel 2004;83 (9):1183–90. [23] Zhang X, Chen Q, Bradford R, Sharifi V, Swithenbank J. Experimental investigation and mathematical modeling of wood combustion in a moving grate boiler. Fuel Process Technol 2010;91(11):1491–9. [24] Nussbaumer T, Kiener M, Horat P. Fluid dynamic optimization of grate boilers with scaled model flow experiments, CFD modeling, and measurements in practice. Biomass Bioenergy 2015;76:11–23. [25] Merrick D. Mathematical models of the thermal decomposition of coal 1. The evolution of volatile matter. Fuel 1983;62:534–9. [26] Merrick D. Mathematical models of the thermal decomposition of coal, 2. Specific heats and heats of reaction. Fuel 1983;62:540–6. [27] Merrick D. Mathematical models of the thermal decomposition of coal: 3. Density, porosity and contraction behavior. Fuel 1983;62:547–52. [28] Atkinson B, Merrick D. Mathematical models of the thermal decomposition of coal. 4. Heat transfer and temperature profiles in coke-oven charge. Fuel 1983;62:553–61. [29] Voller VR, Cross M, Merrick D. Mathematical models of the thermal decomposition of coal. 5. Distribution of gas flow in coke oven charge. Fuel 1983;62:562–6. [30] Klose W, Nowack G. Zum zweidimensionälen Wärmetransport bei chemischkinetisch veränderlichen Stoffwertefunktionen. Wärme- und Stoffübertragung 1986;20:347–60. [31] Patankar VS. Numerical heat transfer and fluid flow. New York: Hemisphere; 1980. [32] Buczynski R, Weber R, Kim R, Schwöppe P. One-dimensional model of heatrecovery, non-recovery coke ovens. Part III: Upper-oven, down-comers and sole-flues. Fuel; 2016 [submitted for publicaion]. [33] Tomeczek J. Zgazowanie we˛gla. Gliwice. Wydawnictwo Politechniki Śla˛skiej; 1991 [in Polish]. [34] Kim R. Untersuchungen zu den Grundsätzen der Wärmetechnik im Heat Recovery –Koksofen. PhD thesis. Clausthal University of Technology. Clausthal-Zellerfeld; 2010. [35] Chiew Y, Glandt E. Simultaneous conduction and radiation in porous and composite materials: effective thermal conductivity. Ind Eng Chem Fundam 1983;22(3):276–82. [36] Hütter ES, Kömle NI. Determination of the radiative contribution to the effective thermal conductivity of a granular medium under vacuum conditions. In: 5th European thermal-sciences conference. Netherlands; 2008. [37] Schotte W. Thermal conductivity of packed beds. AIJChE J 1960(6):63–7. [38] Eisermann W, Johnson P, Cogner WL. Estimating thermodynamic properties of coal, char, tar and ash. Fuel Process Technol 1980(3):39–53. [39] Agroskin A, Goncharov E. The specific heat of coal. Coke Chem. U.S.S.R. 1965:9–14. [40] Agroskin A, Gonczarow L. Thermal capacity and heat of pyrolysis of donbass coal. Koks i Chimija 1970(5):8–13. [41] Gladkov L, Lebedev A. Specific heat of solid fuel and coal dust. Izv VTI 1948 (8):18–20. [42] Kirov N. Specific heats and total heat contents of coals and related materials at elevated temperatures. BCURA Month Bull 1965(29):33–57. [43] Postrzednik S. Analysis of factors determining the course of solid fuel devolatilization and energy consumption. Coke Making Int 1991;3(1):43–7. [44] Wicke M, Peters W. Specific heat capacity. Heat and temperature transfer coefficient of solid fuels. Brennstoffchemie 1968(49):97–102 [in German]. [45] Rycombel S, Długosz A. Appreciation of the thermal effects in the coking process of the selected types of coals. Metalurgia i Odlewnictwo 1979 (5):37–50. [46] Foxwell G. Examination of plasticized coal. J Inst Fuel 1939(12):281–92. [47] Marzec A, Alvarez R, Casal DM, Schulten HR. Basic phenomena responsible for generation of coking pressure: field ionization mass spectrometry studies. Energy Fuels 1995(9):834–40. [48] Walker A. The mechanism of the development of coking pressure. In: Proceedings of 9th international conference on coal science, vol. II. Essen; 1997. p. 601–4. [49] Wasiliev J, Kuzniczenko W. Sposob izmierienia dawljenia razpirania koksujemoj zagruzki w pałupromyszliennoj pieczi. Koks i Chimija 1999 (6):16–21. [50] Barriocanal C, Hays D, Patrick J, Walker A. A laboratory study of the mechanism of coking pressure generation. Fuel 1998(77):729–33. [51] Koch A, Gruber R, Cagniant D, Pajak J, Krzton A, Duchene J. A physicochemical study of carbonization phases. Part I. Tars migration and coking pressure. Fuel Process Technol 1995(45):135–53. [52] Huhn F, Opdenwinkel H, Stoppa H, Strelow F, Strunk J. Influence of raw material and operating parameters onto the internal gas pressure within high capacity coke ovens. In: Extended abstracts and programme of eurocarbon conference – science and technology of carbon, vol. II. Strasbourg; 1998. p. 557–8. [53] Arendt P. Characterization of the swelling behaviour of coals in industrial coke ovens. In: Proceedings of 9th international conference on coal science, vol. II. Essen; 1997. p. 797–800. [54] Heckman H, Klose W. Beeinflussung des Durchgasungsverhaltens der plastischen Zone und des Halbkokses. Gluckauf Forschunghefte 1994 (55):5–14. [55] Sakawa M, Shiraishi K. In-situ observation of coal carbonization process and estimation of coke quality in width direction of oven chamber by computerized tomography (CT). In: 1st international cokemaking conference, vol. I, D2. Essen; 1997. [56] Rozwadowski A. Ciśnienie koksowania mieszanek we˛glowych _ ˛ cych sie˛ właściwościami przygotowanych z komponentów wyraźnie róznia reologicznymi masy plastycznej podczas pirolizy. Gospodarka Surowcami Mineralnymi 2007(21):49–63 [in Polish]. Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086 1029 1030 1031 1032 1033 1034 1035 1036 1037 1038 1039 1040 1041 1042 1043 1044 1045 1046 1047 1048 1049 1050 1051 1052 1053 1054 1055 1056 1057 1058 1059 1060 1061 1062 1063 1064 1065 1066 1067 1068 1069 1070 1071 1072 1073 1074 1075 1076 1077 1078 1079 1080 1081 1082 1083 1084 1085 1086 1087 1088 1089 1090 1091 1092 1093 1094 1095 1096 1097 1098 1099 1100 1101 1102 1103 1104 1105 1106 1107 1108 1109 1110 1111 JFUE 10067 No. of Pages 17, Model 5G 5 February 2016 R. Buczynski et al. / Fuel xxx (2016) xxx–xxx 1112 1113 1114 1115 1116 1117 1118 1119 1120 1121 1122 1123 1124 1125 [57] Rozwadowski A, Strugała A. Badania ciśnienia generowanego w czasie pirolizy przez warstwe˛ plastyczna˛ we˛gli o róźnych właściwościach koksotwórczych. Gospodarka Surowcami Mineralnymi 2006(22):73–81 [in Polish]. [58] Khan M, Gransden J, Price J. Extending battery life by improving coke contraction. In: Proc 1st int mtg lronmaking, Belo Horizonte, Brazil; 2001. p. 239–42. _ [59] Karcz A, Winnicka G. Ciśnienie rozpre˛zania. Cz. II. Wpływ czynników _ technologicznych na ciśnienie rozpre˛zania. Karbo 2001(6):204–8 [in Polish]. [60] Russel CC, Perch M, Smith HB. Measurements of gas pressure within the coal charge in coke ovens. In: Proc blast furn conf AIME, vol. 12; 1953. p. 197–8. [61] Soth GC, Russel CC. Sources of pressure during carbonization of coal. Trans AIME Coal Div 1944(157):281. [62] Nomura S, Thomas M. The effect of swelling pressure during coal carbonization on coke porosity. Fuel 1996;75(2):187–94. 17 [63] Rozwadowski A. The coking pressure of coal blends prepared from components differing markedly in the rheological properties of the plastic mass during pyrolysis. Gospodarka Surowcami Mineralnymi 2007(23):49–63 [in Polish]. [64] Technical coal research. Basic mechanisms of coking pressure generation and the identification of potentially-dangerous coking coals. Luxembourg: European Communities; 1993. [65] Gagarin S. Pressure generation in coking furnaces: a review. Coke Chem 2011;54(6):193–203. [66] British Coal Utilization Research Association. Further studies of coking pressure. Technical Paper No. 5. London; 1952. [67] Buczynski R, Weber R, Kim R, Schwöppe P. One-dimensional model of heatrecovery, non-recovery coke ovens. Part IV: Numerical simulations of the industrial plant. Fuel; 2016 [submitted for publicaion]. 1126 1127 1128 1129 1130 1131 1132 1133 1134 1135 1136 1137 1138 1139 1140 Please cite this article in press as: Buczynski R et al. One-dimensional model of heat-recovery, non-recovery coke ovens. Part II: Coking-bed sub-model. Fuel (2016), http://dx.doi.org/10.1016/j.fuel.2016.01.086