Available online at www.sciencedirect.com
Fuel 87 (2008) 482–490 www.fuelfirst.com
Prediction of ignition behavior in a tangentially fired pulverized coal boiler using CFD T. Asotani b
a,c,*
, T. Yamashita b, H. Tominaga b, Y. Uesugi c, Y. Itaya c, S. Mori
d,1
a Technical Department, Idemitsu Engineering Co., Ltd., 3-1 Nakase, Mihama-ku, Chiba City, Chiba 261-8501, Japan Coal and Environment Research Laboratory, Idemitsu Kosan Co., Ltd., 3-1 Nakasode, Sodegaura, Chiba 299-0267, Japan c Department of Chemical Engineering, Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi 464-8603, Japan d Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi 464-8603, Japan
Received 27 November 2006; received in revised form 9 April 2007; accepted 13 April 2007 Available online 24 May 2007
Abstract Prediction of pulverized coal ignition behavior in a 40 MW tangentially fired commercial boiler is studied. Pulverized coal combustion simulation is performed considering radiation properties of particles. Coal devolatilization and char combustion are modeled and the first order spherical harmonic approximation is used to model the radiative transfer equation. To confirm the accuracy of the simulation method, the results are confirmed by available operating data, design data, and the ignition image in the boiler whose inside is observed by the developed high temperature resistant CCD video camera system. The work indicates that the simulation method can be applied to commercial boilers and predict the ignition behavior with considering not only coal properties but also boiler operating conditions. 2007 Elsevier Ltd. All rights reserved. Keywords: Ignition; Pulverized coal; Tangentially fired boiler; Combustion simulation
1. Introduction Pulverized coal combustion phenomenon can be divided into two steps; devolatilization and char combustion. Ignition, which is a crucial factor for flame stability, the emission of pollutants, boiler design and operating conditions, is a subsequent step of the devolatilization. As the ignition characteristic depends on the type of coal, it should be reflected in boiler design and the operating conditions should also be determined by preliminarily consideration of the compatibility with the design coal. However it is difficult to predict the ignition character of coal only from the proximate and ultimate analysis data. So several experimental methods have been investigated to provide a better *
Corresponding author. Address: Idemitsu Engineering Co., Ltd., 3-1 Nakase, Mihama-ku, Chiba City, Chiba 261-8501, Japan. Tel.: +81 43 296 6951; fax: +81 43 296 6949. E-mail address:
[email protected] (T. Asotani). 1 Emeritus Professor. 0016-2361/$ - see front matter 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.fuel.2007.04.018
evaluation. Although captive particle experiment in TG [1], particle suspension experiment in drop tube furnace [2] and other experiments in a particle suspending apparatus [3] are general experiment methods, these gas flow patterns are not turbulent but laminar and the combustion air temperatures are higher than commercial boiler operating conditions. These differences in combustion conditions between experiments and commercial boilers make it difficult to directly apply the experimental results to commercial boilers, in addition the experimental cost is high. On the other hand, over 20 years CFD (computational fluid dynamics) has been applied to pulverized coal combustion in bench-, pilot-, and commercial-scale furnaces to predict the combustion phenomena in furnaces and substituted for experiment [4–8]. However, the number of applications of CFD to ignition phenomena of pulverized coal combustion is very few. In our previous investigation, an evaluation method of ignition properties by CFD was given and applied to an experimental coal fired furnace (coal feed: 6 kg/h). The results showed the importance of
T. Asotani et al. / Fuel 87 (2008) 482–490
483
Nomenclature a Ac Ap Av C Ec Ev f h Hreac I Ib kv mp n P 2O P0 rp R
surface area of particle (m2/kg) frequency factor of char combustion (kg/ (m2 s Pa)) particle projected area (m2) frequency factor of devolatilization (1/s) weight of char (kg) apparent activation energy of char combustion (kJ/mol) apparent activation energy of devolatilization (kJ/mol) asymmetric factor of particle scattering (–) heat transfer coefficient (W m2) heat generated by char combustion (W) radiation intensity (W m2 sr1) radiation intensity from black body (W m2 sr1) reaction rate of devolatilization (1/s) weight (kg) order of reaction (–) partial pressure of oxygen (Pa) total pressure (Pa) particle reflectivity (=1 ep) (–) universal gas constant (kJ mol1 K1)
the radiation model and the radiation property of coal because the pulverized coal particle temperature is rapidly preheated prior to ignition, mainly by radiation dominating the heat transfer from flames [9]. In this study, to predict the ignition condition in commercial pulverized coal boiler directly, the CFD evaluation method was performed for the commercial 40 MW tangentially fired pulverized coal boiler. The validity of the simulation result was confirmed by available operating data and design data as well as by the observed ignition images in the boiler. 2. Simulation models 2.1. Overview of simulation models The simulations in the present work were done using the commercial CFD software ‘‘FLUENT’’. The standard k–e model and Lagrangian particle tracking method with the random walk model were used for a gas-particle turbulence flow simulation. Although volatile matter contains many kinds of gas, the released volatile matter was represented by a single virtual material. Because the reaction rate of the volatile matter was limited by the turbulent mixing rate of the evolved gas and oxidizer at the near-burner condition, the eddy break-up model was employed for the volatile matter reaction with the following reaction paths. Volatile matter þ a1 O2 ! a2 CO þ a3 H2 O ð1Þ ð2Þ CO þ 0:5O2 ! CO2
s t Ti Uc V V* ep h hR j rs r /
distance along the direction of radiation propagation (m) time (s) temperature (K) ratio of unburned carbon (–) weight of volatile matter (kg) ultimate weight of volatile matter (kg) particle emissivity (–) scattering angle (rad) radiation temperature (=I/4r) (K) absorption coefficient (m1) scattering coefficient (m1) Stefan–Boltzman constant (W m2 K4) circumferential angle (rad)
Subscripts c char d diffusion g gas n number of particle trajectory p particle w wall
The pulverized coal combustion model in this CFD software has been modified as described below. The suitability of the modification was proven by experimental furnace data [10]. 2.2. Coal combustion model The pulverized coal combustion phenomenon should be divided into two processes, namely, devolatilization and char combustion. In FLUENT, coal combustion model consists of the two competing rates devolatilization model [11] and kinetics/diffusion char combustion model [12,13]. The devolatilization model was developed from experiments in an inert atmosphere; on the other hand, the char combustion model was investigated under combustion conditions. Each process was modeled without consideration of the interaction with each other, so the combination requires sequential execution. It should be pointed out, however, that volatile matter evolution and char combustion occurred simultaneously in the initial stage of coal particle combustion [14], and the particle temperature of combusting coal particles is more than 500–600 K higher than the surrounding gas temperature [15], which brings increasing volatile matter evolution. Consequently the sequential execution causes ignition delay in the initial stage of coal particle combustion simulation [10]. In this study, the weight loss of coal particle by devolatilization and char combustion was estimated by global weight loss model proposed by Saito et al. [16] given by
484
T. Asotani et al. / Fuel 87 (2008) 482–490
Eq. (7), and the kinetics/diffusion char combustion model [17] was adopted. The devolatilization rate was represented by a one-step Arrhenius rate law from Eqs. (3) and (4), the char combustion rate was calculated by the difference between the total weight loss and the devolatilization rate. dV ¼ k v ðV V Þ dt k v ¼ Av expðEv =RT p Þ
ð3Þ ð4Þ
The char surface combustion reactions were defined by following reaction paths: Char þ O2 ! CO
ð5Þ
CO þ 0:5O2 ! CO2 n dC 6 d p32 Ec P O2 ¼ Ac exp C dt c qp RT p P0 0:75 dC ½ðT p þ T g Þ=2 ¼ d 2p P O2 2:53e 12 dt d dp dC 1 ¼ 1 dt þ dC1 ðdC ð dt Þd dt Þc
ð6Þ ð7Þ ð8Þ ð9Þ
The devolatilization and char combustion were calculated simultaneously during pulverized coal particle trajectory calculations. These devolatilization and char combustion parameters were measured in an electrically-heated drop tube furnace, the former was measured under inert atmosphere conditions, the latter was measured under combustion conditions. The coal properties and the parameters of devolatilization and char combustion are listed in Table 1.
Table 1 Coal properties and parameters Proximate analysis (wt%, ad) Moisture Ash Volatile matters Fixed carbon
3.4 14.3 31.6 50.7
Ultimate analysis (wt%, daf) Carbon Hydrogen Nitrogen Sulfur Oxygen
83.48 5.46 2.00 0.43 8.63
Calorific value (MJ/kg, ad)
28,424
Devolatilization parameters Av (l/s) Ev (kJ/mol) V* (–)
2021 31.1 0.43
Char combustion parameters Ac (kg/m2/s) Ec (kJ/mol) n (–)
10.2 57.1 0.5
Particle diameter distribution Under #200 (wt%) Rosin–Rammler distribution constant (–)
80 1.19
The particle diameters are fitted to the Rosin–Rammler distribution, with representative diameters 4, 11, 17, 24, 32, 41, 52, 65, 85, and 125 lm. 2.3. Radiation model Considering radiation intensity along a path s through an absorbing, emitting and scattering mixture of coal combustion gases and particles in local thermodynamics equilibrium, the radiative transfer equation (RTE) can be written as Z dI rs ¼ ðj þ rs ÞI þ jI b þ Ið/ÞpðhÞd/ ð10Þ ds 4p 4p where I(/) is the incident radiative intensity from angle / and p(h) is the particle scattering phase function. In highly forward scattering media, the effect of scattering on radiative transfer can be neglected, and the results of radiative heat flux obtained by spherical harmonics (PN) are in good agreement with numerically exact results [18]. The authors examined the directional behavior of scattering in coal ash particle clouds by using FT-IR. The results indicated that the contribution of scattering by coal particles can be ignored, and the forward scattering dominates the particle scattering in the radiative heat transfer in the cloud [19]. The PN approximation can be applied to the modeling of the RTE, and can be used with finite-difference schemes required for flow-field calculations and yield reasonably accurate radiative heat flux predictions, especially for optical thickness greater than unity [20]. In this study, the RTE is modeled using the first order PN approximation for the angular-dependency of the radiative flux. In this approximation, the equation for the radiation temperature, hR is i Xh 1 1 4 r hR j h4R T 4 ep Apn ðh4R T 4pn Þ ¼ 0 ð11Þ 3 jR n X jR ¼ j þ Apn ep þ rp ð1 f Þ ð12Þ n
The highly forward scattering of radiation by pulverized coal particles is accounted by delta-Eddington phase function approximation [21], which is given by asymmetric factor: f in Eq. (12) and set f = 1.0. The absorption coefficient of the gas phase is determined by the weighted-sum-ofgray-gases model (WSGGM) [22]. The soot concentration is assumed to be uniform and constant (103 kg/m3) in the furnace, and the absorption coefficient of gas phase is corrected by the Taylor–Foster model [23]. In the mixture of coal combustion gases and particles, radiative transfer is mainly caused by particle clouds, so the treatment of particle emissivity ep is essential in the pulverized coal combustion model. In the pulverized coal experimental furnace (coal feed 6 kg/h), the parametric study indicates that particle emissivity influences the temperature distribution especially near the burner. With increasing of particle emissivity, the difference between the calculated and measured temperature at the near-burner measurement point (about
T. Asotani et al. / Fuel 87 (2008) 482–490
485
130 mm from burner nozzle) decreases from 100 K at ep = 0.4 to 10 K at ep = 1.0. The calculated temperature is in agreement with the measured one when ep is assumed to be a function of unburned carbon rate [9]. In this study, the calculations included the effects of particle radiation; the particle emissivity ep was assumed to be a function of unburned carbon in particle Uc [24], i.e. ep ¼ 0:4 U c þ 0:6
ð13Þ
During combustion, the particle temperature is estimated by the following heat balance equation, and application of the calculation method presented by Mitchell [25]. dmp H reac þ Ap ep r T 4p h4R qp ¼ hAp ðT p T g Þ ð14Þ dt 3. Simulation and observation methods 3.1. Simulation conditions The simulation model described above was applied to a 40 MW commercial pulverized coal fired boiler, the furnace has dimensions 5.5 · 5.5 · 17.5 m, with tangential fired burners located in the corners. The furnace geometry and burner arrangement are shown in Fig. 1, and the finite-difference grid is depicted in Fig. 2. There are 54 · 40 cells in a cross-section, totally 54 · 40 · 129 = 278,640 hexahedral cells. This boiler has two types of coal burner, one is a fuel rich burner (CONC) the other is a fuel lean burner (WEAK). The geometry of these burners is meshed in detail in the analysis. The burner size, the boundary condi-
Fig. 2. Grid formation at the cross-sections, furnace center and CONC burner level.
tions and the operating conditions are listed in Tables 2 and 3. Divided 10 particle diameters of the pulverized coal are initiated at 10 injection points in each coal burner, so that a total of 1200 particles are calculated. 3.2. Observation method In general, it is difficult to directly observe the ignition condition especially in tangentially fired boilers because of the constructional limitations. In order to observe the ignition behavior clearly, the high temperature resistant CCD video camera system was developed [26]. The system consists of three concentric water-cooled cylinders, at the tip of which a high performance CCD and a ceramic mirror Table 2 Boundary conditions
Fig. 1. Boiler and burner geometry (CONC, WEAK: fuel rich and lean coal burner; AUX, OIL, SGR: combustion air feed port; OFA: over fire air feed port; U: upper; L: lower).
Burner air port
Width (mm)
Height (mm)
Velocity (m/s)
OFA AUX OIL SGR CONC WEAK
250 270 300 250 400 300
80 40 240 40 250 150
22.4 6.1 6.2 9.0 8.6 10.3
Radiation boundary conditions Tw (K) 700 ew (–) 0.6
486
T. Asotani et al. / Fuel 87 (2008) 482–490
Table 3 Boiler operating conditions
Table 4 Validation of simulation results
Heat value (MW) Total coal flow (kg/s) Excess air (%) Combustion air temperature (K)
40 3.67 18 622
Coal burner conditions Air and coal temperature (K) Air/coal (kg/kg) CONC burner WEAK burner Velocity (m/s) CONC burner WEAK burner
350
Burner
Unburned carbon in ash (wt%)
Upper
CONC WEAK
15.03 6.47
2.10 1.86
Lower
CONC WEAK
6.00 0.93
4.05 8.18
9.82 10.30
4.15 –
0.89 5.05
Total (calculated) Measured
17.7 25.5
Outlet temperature (K)
are mounted. Purging air that prevents exposure to hot flue gas and ash, and capture of the image through the ceramic mirror enables direct observation in the furnace. The high temperature resistant CCD video camera was set up at the peephole located under the lower burner as depicted in Fig. 1, and the view angle is showed by the dotted line. 4. Results and discussion
Particle residence time (s)
Calculated result 1350
Boiler design 1373
devolatilization and char combustion at particle surface) are executed alternately until the total mass and heat are balanced. The simulated result is validated by the operating data and the boiler design parameter, including the average gas temperature and the unburned carbon in ash at the furnace outlet. As shown in Table 4, good agreement with the measured and the design data confirms the accuracy of the simulation method.
4.1. Validation of simulation results The calculation starts by solving only the isothermal gas flow without particles. After the convergence, the particle trajectories, the chemical reactions and the enthalpy are solved, finally the radiation transfer calculation is included. The gas phase calculation (flow, heat, mass and reaction) and the pulverized coal particle calculation (trajectory,
4.2. Gas flow and temperature distribution and particle trajectory The gas velocity and the temperature distribution at the cross-sections between the opposed burners and the crosssections along the furnace height are shown in Fig. 3. In this figure, a tangentially-swirling flow, which is created
Fig. 3. Velocity vector and temperature distribution in the furnace.
T. Asotani et al. / Fuel 87 (2008) 482–490
by the air injection from the burner port, produces a socalled ‘‘fire ball’’ in the furnace center. The temperature distributions show the flame heat release variation with furnace height as a change of gas temperature, the maximum temperature of about 1850 K is observed at the height just above the upper-WEAK burner ports. Heat transfer to the furnace wall causes a temperature decrease in the flue gas flow, consequently the average temperature at the outlet becomes 1350 K. Fig. 4 shows mean particle trajectories (shaded by particle residence time) fed from each burner without turbulence dispersion. It is clear that many particles from the lower burner circulate and stay around the bottom of the furnace and eventually pass through the furnace center. On the other hand, particles from the upper burner could not
487
reach the furnace center and instead flow around the surface of ‘‘fire ball’’, as a result, the particle residence time is short and the unburned carbon in ash is relatively high compared with particles from the lower burner, as shown in Table 4. Comparing the particles fed from the WEAK and CONC burners, particles from WEAK burner go forward to the furnace center because the injected velocity is higher than that of CONC burner. 4.3. Ignition of pulverized coal clouds Near-burner images captured by the CCD video camera system are shown in Fig. 5. In the image of the CONC burner, ignition begins from the right side of the picture that faces the furnace center and the ignition boundary
Fig. 4. Mean particle trajectories fed from each burner without turbulence dispersion (shaded by particle residence time).
Fig. 5. Near-burner images captured by the CCD video camera system.
488
T. Asotani et al. / Fuel 87 (2008) 482–490
approaches the burner face. On the other hand, ignition does not occur near the WEAK burner. Simulated distribu-
tions of volatile matter, CO concentration and temperature are shown in Figs. 6–8. Near the CONC burner, the
Fig. 6. Calculated distributions of volatile matter concentration (for comparison, above 5% is displayed in same color).
Fig. 7. Calculated distributions of CO concentration (for comparison, above 5% is displayed in same color).
Fig. 8. Calculated distributions of temperature.
T. Asotani et al. / Fuel 87 (2008) 482–490
489
Fig. 9. Calculated distributions of gas velocity.
evolved volatile matter reacts with oxygen in gas phase and rapid temperature rise occurs, in particular the distribution of CO concentration looks like the shape of ignition boundary in the video camera image. On the other hand, simulation results for the WEAK burner show low volatile matter and CO concentrations, and no temperature rise caused by ignition, similarly to the observed images. To illustrate the differences between the CONC and the WEAK burner, simulated distributions of gas velocity and coal particle concentration are shown in Figs. 9 and 10. As a result of the boundary conditions of velocity and air/coal ratio shown in Table 3, near the WEAK burner the gas velocity is higher and the pulverized coal particle concentration is lower than those of the CONC burner. These are the cause of the ignition condition difference. In general, pulverized coal clouds fed from burner receive the radiative energy from flames in furnace center, which causes rapid increasing of coal particle temperature with devolatilization. In Fig. 11, radiation temperature, which indicates the intensity of radiation energy, are compared for each burner. The radiation temperature near the
CONC burner is lower than that of the WEAK burner. The difference indicates that the radiation energy that is propagated from the furnace center is absorbed by denser pulverized coal clouds before it reaches the CONC burner, which causes rapid increasing of pulverized coal particle temperature and hence devolatilization. Furthermore, the high concentration of the evolved volatile matter (because of the low gas velocity) results in reaction with oxygen in gas phase. In these simulation results, pulverized coals fed from CONC burner receive radiative energy which is propagated from the furnace center, the devolatilization occurs in the near-burner region and the volatile matter distribution looks like the shape of ignition boundary in the video camera image. The simulated temperature distribution shows the rapid temperature rise in this area, which indicates a stable flame, the same as observed by the CCD video camera system. On the other hand, simulation results for the WEAK burner show low volatile matter concentration and no temperature rise caused by ignition, similar to the observed images.
Fig. 10. Calculated distributions of pulverized coal particle concentration.
490
T. Asotani et al. / Fuel 87 (2008) 482–490
Fig. 11. Calculated distributions of radiation temperature (=I/4r).
5. Conclusions Simulations were conducted to predict the ignition condition of pulverized coal clouds in a 40 MW commercial tangentially fired boiler. Good agreement of simulated results with the observed images, the operation and design data indicates the accuracy of the simulation method and the application to this commercial boiler. It is proven that this simulation method can evaluate pulverized coal ignition behavior directly for commercial boilers by considering not only coal properties but also boiler operating conditions, and the substitution for experiment will bring reduction of the cost. References [1] Tognotti L, Malotti A, Petarca L, Zanelli S. Measurement of ignition temperature of coal particles using a thermogravimetric technique. Combust Sci Tech 1985;44:15–28. [2] Wall TF, Phong-anant D, Gururajan VS, Wibberley LJ, Tate A, Lucas J. Indicators of ignition for clouds of pulverized coal. Combust Flame 1988;72:111–8. [3] McLean WJ, Hardesty DR, Pohl JH. Direct observations of devolatilizating pulverized coal particles in a combustion environment. In: 18th Symposium (international) on combustion, 1981. p. 1239–48. [4] Fiveland WA, Wessel RA. Numerical model for predicting performance of three dimensional pulverize-fuel fired furnaces. J Eng Gas Turb Power 1988;110:117–26. [5] Lockwood FC, Mahmud T. The prediction of swirl burner pulverized coal flames. In: 22nd Symposium (international) on combustion, Pittsburgh, USA, 1988. p. 165–73. [6] Boyd RK, Kent JH. Three-dimensional furnace computer modelling. In: 21st Symposium (international) on combustion, Munich, West Germany, 1986. p. 265–74. [7] Tominaga H, Sato M, Kambara S, Inouch K, Yamada T. Simulation of pulverized coal combustion: Prediction of unburned carbon at furnace outlet. In: Proceedings of the 7th international conference on coal science, Newcastle, UK, 1991. p. 301–10. [8] Smith JS, Smith PJ, Hill SC. Parametric sensitivity study of a CFDbased coal combustion model. Am Inst Chem Eng J 1993;39:1668–79. [9] Asotani T, Yamashita T, Tominaga H, Uesugi Y, Itaya Y, Mori S. Effect of radiative heat transfer properties on combustion simulation for pulverized coal experimental furnace. In: 13th International heat transfer conference, Sydney, NSW Australia, 2006; COM-15.
[10] Tominaga H, Harada M, Ando T, Suzuki Y. The development of pulverized coal combustion simulator. J NIRE 1997;6:333–43. [11] Kobayashi H, Howard JB, Sarofim AF. Coal devolatilization at high temperatures. In: 16th Symposium (international) on combustion, 1977. p. 411–25. [12] Baum MM, Street PJ. Predicting the combustion behavior of coal particles. Combust Sci Tech 1971:231–43. [13] Field MA. Rate of combustion of size-graded fractions of char from a low rank coal between 1200–2000 K. Combust Flame 1969:237–52. [14] Howard JB, Essnhigh RH. Mechanism of solid-particle combustion with simultaneous gas-phase volatiles combustion. In: 11th Symposium (international) on combustion, 1982. p. 399. [15] Timothy LD, Sarofirm AF. Characteristics of single particle coal combustion. In: 19th symposium (international) on combustion, 1982. p. 1123. [16] Saito M, Sadakata M, Sakai T. Measurements of surface combustion rate of single coal particles in laminar flow furnace. Combust Sci Technol 1987;51:109–28. [17] Mulcahy MFR, Smith IW. Kinetics of combustion of pulverized fuels: a review of theory and experiment. Rev Pure Appl Chem 1969;19:81–108. [18] Liu F, Garbett ES, Swithebank J. Effect of an anisotropic scattering on radiative transfer using P1-approximation. Int J Mass Transfer 1992;35:2491–9. [19] Itaya Y, Nishio N, Hatano S, Kobayshi N, Kobayashi J, Mori S. Near-infrared radiation and scattering properties of coal fly ash particles cloud. J Chem Eng Jpn 2006;39:718–23. [20] Menguc MP, Viskanta R. Radiative transfer in three-dimensional rectangular enclosures containing inhomogeneous, anisotropically scattering media. J Quanti Spectrosc Radiat Transfer 1985;33:533–49. [21] Joseph JH, Wiscombe WJ, Weinman J. The delta-Eddington approximation for radiative heat transfer. J Atomos Sci 1976;33: 2452–9. [22] Smith TF, Shen ZF, Friedman JN. Evaluation of coefficients for the weighted sum of gray gases model. J Heat Transfer 1982;104:602–8. [23] Taylor PB, Foster PJ. Some gray gas weighting coefficients for CO2H2O-soot mixtures. Int J Heat Mass Transfer 1974;18:1331–2. [24] Lockwood FC, Rizvi SM, Shah NG. Comparative predictive experience of coal firing. Proc Inst Mech Eng 1986;200:79–87. [25] Mitchell RE. On the temperature and reaction rate of burning pulverized fuels. In: 19th symposium (international) on combustion 1982. p. 1113–22. [26] Asotani T, Yamashita T, Tachibana K. Direct observation of ash deposition in pulverized coal fired boiler using CCD camera and model development for slagging prediction. In: 6th World conference on experimental heat transfer, fluid mechanics, and thermodynamics, Matsushima Miyagi Japan, 2005, 9-a-15.