Next Article in Journal
Critical Clearing Time and Wind Power in Small Isolated Power Systems Considering Inertia Emulation
Next Article in Special Issue
On the Heat Flux Vector and Thermal Conductivity of Slags: A Brief Review
Previous Article in Journal
The Influence of Non-Uniform High Heat Flux on Thermal Stress of Thermoelectric Power Generator
Previous Article in Special Issue
Effect of Heterogeneity in Coal Ash Chemical Composition on the Onset of Conditions Favorable for Agglomeration in Fluid Beds
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Review

Modelling Underground Coal Gasification—A Review

1
Department of Chemical and Materials Engineering, University of Alberta, Donadeo Innovation Centre for Engineering, 9211–116 Street NW, Edmonton, AB T6G 1H9, Canada
2
School of Mining and Petroleum Engineering, University of Alberta, Edmonton, T6G 1H9, AB, Canada
*
Author to whom correspondence should be addressed.
Energies 2015, 8(11), 12603-12668; https://doi.org/10.3390/en81112331
Submission received: 18 August 2015 / Revised: 25 October 2015 / Accepted: 27 October 2015 / Published: 6 November 2015
(This article belongs to the Special Issue Recent Advances in Coal Combustion and Gasification)

Abstract

:
The technical feasibility of underground coal gasification (UCG) has been established through many field trials and laboratory-scale experiments over the past decades. However, the UCG is site specific and the commercialization of UCG is being hindered due to the lack of complete information for a specific site of operation. Since conducting UCG trials and data extraction are costly and difficult, modeling has been an important part of UCG study to predict the effect of various physical and operating parameters on the performance of the process. Over the years, various models have been developed in order to improve the understanding of the UCG process. This article reviews the approaches, key concepts, assumptions, and limitations of various forward gasification UCG models for cavity growth and product gas recovery. However, emphasis is given to the most important models, such as packed bed models, the channel model, and the coal slab model. In addition, because of the integral part of the main models, various sub-models such as drying and pyrolysis are also included in this review. The aim of this study is to provide an overview of the various simulation methodologies and sub-models in order to enhance the understanding of the critical aspects of the UCG process.

1. Introduction

Coal is the most abundant and widespread type of fossil fuel, with reserves of economically recoverable coal estimated to be about 900 gigatons (GT), 64% of major fossil fuels [1]. However, the total amount of coal potential including the economically non-recoverable part is reported to be 18 terratons [2], most of which are too deep or too costly to be mined by traditional methods. With the current trends of energy consumption, and the rise of developing countries, the total fossil fuel reserve is expected to be the same while the total energy consumption is increasing around the world. In order to increase the recoverable fossil fuel reserve, it is important to find appropriate technology to extract the non-recoverable part of the reserves economically. Underground coal gasification is a method of in-situ conversion of deep un-minable coal into gaseous products with a higher heating value through controlled combustion with oxygen/air and gasification with steam. The UCG product gas, called synthesis gas (syngas), is mainly a mixture of CO, H2, CH4, and CO2 with trace amounts of ethane and other components. The development of the underground (in-situ) coal gasification (UCG) technique is paving the way to access deep coal seams. Early studies suggest that the use of UCG could potentially increase the world’s coal reserves by as much as 600 GT [3], which represents a 70% increase [2].

1.1. Environmental and Safety Features of UCG and Coal Mining

The technology of coal mining and the use of coal have been affected by the pollution caused by its transport, storage, and combustion [4]. In addition, coal-fired power plants have higher greenhouse gas (GHG) emission per generated energy due to lower efficiency of the process and a higher carbon-to-hydrogen ratio of coal. However, the most alarming issue is the worldwide death fatalities and injuries per year in the coal mining industries. It has been reported that there was a total 36,910 deaths in a seven-year period, from 2000 to 2005, in the coal mining industries of China alone [5]. On the other hand, the USA has reduced their death fatalities after ensuring safe work environments in the current century; however, the death fatalities of the USA have still numbered more than 300 since 2003 [6]. Considering the death fatalities alone, UCG is the safest energy harness process. Moreover, UCG has several advantages over the conventional gasification process. Examples of such benefits include [7,8,9]:
-
Elimination of coal mining and thus avoiding a dangerous practice even with current equipment and safety regulations.
-
Lower capital investment due to the absence of surface gasification units.
-
Elimination of ash handling operations by keeping the ash underground.
-
Lower water consumption.
-
Generation of possible sites for CO2 sequestration.
From the economic and environmental point of view for large-scale gas production, UCG is a highly promising technology [7,10]; however, UCG introduces some challenges that should be addressed before the process can be adopted on a large scale. Some of the main challenges are process stability, aquifer contamination, and subsidence [11]. However, the potential risks of subsidence and aquifer contamination can be avoided through careful site selection and adoption of best management practices for operations [7,8,11,12].

1.2. UCG Process Overview

The UCG process involves four main steps, namely drilling, linking, ignition, and gasification. The basic layout of the process requires one borehole for the injection of gases and one for production; however, over the years three standard configurations of wells have evolved. They are:
-
Linked vertical wells (LVW);
-
Controlled retraction injection point (CRIP);
-
Steeply dipping seams (SDS).
Figure 1 shows a sketch of linked vertical wells arrangements for the UCG process. In the first step, injection and production wells apart from each other are drilled at some specified distance from the surface into the coal seam. In the following step, a permeable channel (link) is created between two wells which allows the sufficient flow of gases in a large-scale operation. It is necessary to ignite the coal initially in order to raise the coal temperature using combustible gas and oxygen mixture with the help of a suitable ignition source. The duration of this combustion stage depends on the seam properties as well as the operating conditions. After this stage, the mixture of air/oxygen and steam are injected into the seam to gasify the coal. The main product of gasification, syngas, is then recovered from the production well. As the coal is being consumed, ash falls down and a void space is created, which is referred to as a cavity. When the cavity reaches the production well, there will be less contact between the injection gases and the virgin coal; thus, the quality of produced gases drops significantly. Therefore, fresh coal should be made available by moving the location of the injection gases for further continuation of the process. After the completion of gasification, steam and/or water are usually flushed to the cavity in order to remove the pollutants from the coal seam to prevent them from diffusing into surrounding aquifers.
Figure 1. Schematic diagram (not to scale) of the UCG process for linked vertical wells (LVW).
Figure 1. Schematic diagram (not to scale) of the UCG process for linked vertical wells (LVW).
Energies 08 12331 g001

1.3. Chemical Reactions in UCG

The UCG process occurs in three different zones: drying, pyrolysis, and combustion and gasification of solid char. In the drying zone, wet coal is heated and converted into dry coal by removing moisture attached to it. Upon further heating of the coal, the pyrolysis reactions begin at temperatures around 350–400 °C [13] and coal loses its weight, generating volatile matters and the solid known as char. Finally, the char reacts with the injected/pyrolyzed gases to produce the syngas. The three essential types of char layer reactions which govern the UCG product gas composition are: oxidation, reduction, and pyrolysis. Many reactions occur during this process, but the most important reactions, which are considered in most of the models, are summarized in Table 1.
Table 1. Main reactions in the gasification process.
Table 1. Main reactions in the gasification process.
ReactionReaction No.RepresentationHeat of Reaction (kJ/mol)
Drying1Coal → Dry Coal+ H2O+40
Pyrolysis2Dry Coal → Char + Volatiles0
Coal combustion3C (Char) + O2 → CO2−393
CO2 gasification4C (Char) + CO2 → 2CO+172
Steam gasification5C(Char) + H2O → CO+ H2+131
Methanation6C(Char) + 2H2 → CH4−75
Water-gas shift7CO + H2O ↔ CO2+ H2−41
Gas phase oxidations8CO + 0.5O2 → CO2−111
9H2 + 0.5O2 → H2O−242
10CH4 + 2O2→ CO2 + 2H2O−802
The reactivity of the char to O2, H2O, CO2, and H2 determines the rates at which the desired products are formed [9]. Reactions 3, 5, and 8 are the main chemical reactions considered for both shallow and deep coal gasification processes. The hydrogasification (reaction 6) is favourable at a high hydrogen pressure. In the UCG, at low pressure, this reaction is not significant. Reactions 3–6 take place on the wall of coal seams, while reactions 7 and 8 occur at the gaseous stage [14]. Different reactions are favourable at different temperature ranges. Figure 2 shows the equilibrium constants as a function of temperature for some of the important reactions in the gasification process. As the equilibrium constants for all the reactions shown in Figure 2 are at unity at approximately 1000 K, the temperature becomes the key factor in determining the shift of chemical equilibria towards the product or reactions by deciding whether it is an endothermic or exothermic reaction. For example, at temperatures above 1000 K, the formations of H2 and CO are favourable since reaction 4 and reaction 5 are endothermic [15]. Similarly, other reactions seem to be favoured at temperatures below 1000 K. Reactions involved in UCG are not fundamentally different from the reactions of coal and char in the surface gasifiers. In gasification reactions, the combustion of coal and volatiles provides the necessary heat for endothermic gasification reactions.
Figure 2. Equilibrium constant Kp as a function of temperature for some reactions in coal gasification [15].
Figure 2. Equilibrium constant Kp as a function of temperature for some reactions in coal gasification [15].
Energies 08 12331 g002

1.4. Brief History of UCG Practice

Although the idea of the UCG was first mentioned by Sir William Siemens of Great Britain in 1868 [16,17,18] and the idea was consolidated (1988–1899) independently by the famous Russian chemist Dmitri I. Mendeleev, who also pointed out the economic benefit of UCG over conventional mining [16], most credit goes to an American chemist, A.G. Betts, due to his detailed technical contents and engineering drawing on the UCG method that closely resemble modern approaches presented in his three UCG patents received during 1909–1910 [16,17]. Later on, in 1912, famous British chemist Sir William Ramsay expanded the ideas of Betts and proposed gasifying coal underground as a way to avoid air pollution and as a solution for the mining worker [16]. His strong advocacy for the development of this technology created an international surge of interest in UCG that eventually helped him to culminate the first ever UCG experiment being carried out at Hett Hill near Tursdale Colliery in County Durham, North East England, in 1912 [16,18]. Unfortunately, the experiments did not proceed due to the outbreak of World War I in 1914 and Ramsay’s death in 1916. A detailed description of early ideas in underground coal gasification and their evolution is given in the excellent historic review of UCG by Klimenko [16].
After successful operations of UCG throughout the 1930s and 1960s, the process was left out mostly because of the low prices of gas and oil. Due to increased energy demands and stringent environmental regulations, UCG has regained its popularity in recent years as an immediate alternative for conventional coal gasification plants. Over 50 pilot UCG plants have been conducted worldwide since the 1930s. Most of these developments have been concentrated in the former USSR, Europe, the USA, and China. The former USSR was the first country to operate UCG successfully and they gained extensive knowledge in the UCG process. However, their knowledge was scarcely available until the 1960s. Most of the information on Russian tests was available in the Russian literature. As a result, most of the field trials conducted by other countries before the 1960s, without the benefit of detailed data on Soviet field tests, did little to advance the process and left an impression that UCG was not a viable process [19]. Later on, many Russian articles on UCG were translated into English. As a result, the volumes of works on UCG carried out by the former USSR were exposed to the rest of the world. After getting a complete picture of the Soviet UCG experience with the large-scale translation of Russian technical literature, Gregg et al. [19] concluded that “the amount of effort expended and degree of success achieved by the Soviets (having constructed and operated several commercial plants) far exceeds the summation of efforts by other nations”. Translations of a Russian book and summaries of the translations of Russian articles and specific experiences on UCG can be found in several references: Skafa [20], Gregg et al. [21], Derbin et al. [22]. After gaining knowledge from the Soviets' field practices, extensive field trials were performed in the USA and large-scale pilots proved the basic feasibility of UCG [23]. Earlier, the field tests were limited to shallow depth, following the Soviet practice, but over time, with the improvements in drilling technology, the depth of UCG operations has increased markedly as can be observed in recent European field trials [11]. The underground preparation of the seam and gasification are performed using drilled holes. Earlier, linked vertical wells (LVW) and steeply dipping seams (SDS) were the most popular drilling configurations for horizontal and steeply dipping seams, respectively. However, recently, the controlled retraction injection point (CRIP) invented in the Centralia series of UCG trials in the 1980s has gained considerable attraction due to its greater control and improved overall efficiency [2,10]. Establishing a highly permeable path between the injection and production holes is considered as the single most important task prior to gasification because a highly permeable narrow channel can control gas leakage, water intrusion, liquid condensation, phenol contamination, and oxidant bypass during UCG operation [19]. Directional drilling, reverse combustion, and hydrofracturing are the popular linking methods currently [2,19]. After the invention of reverse combustion linking (RCL) in 1941 by former Soviet researchers [24], the technique was widely used by many countries in UCG field trials, including the former USSR. Reverse combustion linking (RCL) is the opposite of forward combustion linking (FCL) where the coal is ignited at the base of the injection well and the flame propagates towards the production well, i.e., in the direction opposite to the gas flow within the channel [24,25]. As compared to FCL, RCL is the preferable linking option in terms of channel dimension, consumption of coal per unit of the link, linking speed, operational stability, and controllability [24]. However, FCL can also be used for widening the hydraulic links between the wells which exist naturally or have been preliminarily established using, for example, hydrofracturing [24]. Hydrofracturing is employed to establish a linking channel in coal seam of low permeability where RCL alone is not a suitable option initially due to in-situ coal properties and the fact that directional drilling is too sophisticated for that situation. The effectiveness of hydrofracturing is independent of both the electrical properties and the natural permeability of the coal; however, its lack of control limits its use primarily for steeply dipping beds where it is not essential to ensure that the cracks form at the bottom of the seam or that they form parallel or perpendicular to the bedding planes [19]. Directional drilling is one technique that is relatively independent of the coal properties. The in-seam drilling of coal seams has been part of coal exploitation since at least the 1940s in Soviet plants [19]. Due to continuous and significant advances in the techniques used for directional drilling, this method is being currently considered as the most reliable way of establishing the link between wells.
UCG field tests have been carried in different coal seams worldwide, using different techniques in order to evaluate the effects of various factors of overall UCG performance. Some significant field trials worldwide are summarised in Table 2. A quick look at Table 2 indicates a significant increased heating value of the product gas for an injection gas of oxygen/steam mixture over air. Moreover, the overall performance is largely affected by coal properties (i.e., ranks and chemical composition), coal seam geology (i.e., depth, thickness), and drilling configurations. The dependence of coal properties, coal seam geology, and overburden properties (geology and hydrology) indicates that the performance of UCG is site-specific. The site selection is perhaps the single most important decision that will determine the technical and economic feasibility of coal gasification [15]. However, for same site, seam preparation and operating conditions (i.e., temperature, pressure) are the key factors for the successful performance of UCG.
Table 2. Summary of significant underground coal gasification (UCG) field trials worldwide [2,15,26,27,28,29].
Table 2. Summary of significant underground coal gasification (UCG) field trials worldwide [2,15,26,27,28,29].
CountrySiteStart-Up YearCoal RankWell Config-urationInjected GasSeam Depth (m)Seam Thickness (m)Gas Calorific Value (MJ/m3)
USSR (former)Gorlovka (Ukraine)1937BSDSAir401.93.95
O2/H2O401.95.82–10.34
Podmoskova/Tula (Russia)1947LLVWAir5533.4
Yuzhno-Abinsk1955BLVWAir~10024.1
Angren (Uzbekistan)1961LLVWAir11043.5–3.9
USAHanna II1975BLVWAir8594.2
Hoe Creek I1976BLVWAir4083.6
Hoe Creek IIB1977BLVWO2/H2O4089.0
Price town 11979BLVWAir27026.1
Rawlins IA1979SBSDSAir105185.6
Rawlins IIB1979SBSDSO2/H2O105188.1
Centralia A1984SBCRIPO2/H2O7569.7
Centralia B1984SBLVWO2/H2O7568.4
Rocky Mountain IA1987SBCRIPO2/H2O11079.5
Rocky Mountain IB1987SBLVWO2/H2O11078.8
Alaska SHR2012L/SB Air50–6501–123.3–5
BelgiumThulin1986ALVWAir86067.0
SpainEl Tremedal1997SBCRIPO2/H2O58026.6
AustraliaChinchilla2000SBLVWAir140106.6
Bloodwood2011 CRIPAir2008-106.4
ChinaSuncan2002BTunnelAir8028.5
South A.Majuba2007BLVWAir2853.5–4.56.2
B = Bituminous, SB = Subbituminous, A = Anthracite, L = Lignite.

2. UCG Modeling

Although a number of UCG field trials have been performed, the information on the detailed UCG process is very limited because of the high cost of extracting data as well as the difficulty in controlling the operating variables. In addition, there is a limited scope of adjusting parameters due to the site-specific nature of the performance of UCG. As a result, due to relative economic consideration, several laboratory-scale coal block experiments have been reported [30,31]. Table 3 represents major laboratory studies on UCG. Feed temperature, feed flow rate, steam-to-oxygen ratio, operating pressure, combustion time, total operation time, and distance between wells are varied in order to observe their effects on the product gas composition and flow rate and gas outlet temperature. As well, the temperature distribution in the coal seam during gasification, the final cavity dimension, and the rate of recession of the coal surface due to reaction are also observed and, finally, the optimum operating conditions for the expected UCG performance are determined.
Table 3. Major UCG experiments on coal blocks [27].
Table 3. Major UCG experiments on coal blocks [27].
ResearchersCoal TypeDimensions of Seam (Length × Width × Height)Major Observations
Yeary and Riggs [32]Lignite- Sub-bituminous25 cm × 5 cm × 25 cmLateral cavity growth
Poon [33]Lignite9 cm × 2.5 cmLateral cavity growth
Daggupati et al. [34]Lignite30-38 cm × 20cm × 25cmCavity growth and shape
Yang [35]High-volatile bituminous6.8 m × 0.25 m ×1.1 m (68°) 7.5 m × 2 m × 1 m (65°)Two-phase UCG in inclined abandoned coal mines
Prabu and Jaynati [36]Camphor, Wood, Lignite24 cm × 8–10 cm × 10–12 cmCavity shape
Stanczyk et al. [37]Lignite- Hard Coal2.5m × 0.7m × 0.7mTwo-phase UCG in lignite
The generalizing data from the lab-scale experiments should still be done with extreme care since the whole process may not be properly represented in lab-scale experiments. It may not be possible to represent the entire phenomenon in laboratory-scale experiments [38,39]. The coal seam is a three-dimensional reactor without specific geometry and the flow pattern inside the UCG cavity is highly non-ideal. The complexity increases further because of several other complex physical and chemical processes occurring simultaneously, such as chemical reactions and kinetics, transport phenomena (i.e., heat and mass), water intrusion from surrounding strata, thermo-mechanical processes related to the structure of the seam, and other geological aspects. Both field tests and lab-scale experiments are not sufficient to provide a detailed quantitative description due to the demand of extensive instrumentation and trials. The difficulties encountered in controlling and operating the UCG process serve to be the incentive to develop quantitative models in order to predict the effects of various physical and operating parameters on the performance of the process and to extrapolate the understanding to other operating conditions and coals. Prior to 1975, the development of UCG models was very limited [40], and focused only isolated cases, such as heat transfer, mass transfer, combustion rate, etc. However, after gaining a better understanding of the UCG process, the intensity of developing a global model increased significantly. A full UCG model that describes the complete mechanisms would typically include a number of the following sub-models:
-
Injection/production linkage sub-model.
-
UCG reactor sub-model.
-
Ground water hydrology.
-
Ground subsidence model.
-
Surface facility models.
Combining all these sub-models would theoretically give an exact description of the process; however, building such a model is not a trivial task. Therefore, all previous models have focused on studying these aspects separately through several simplifying assumptions. There are lots of models developed thus far; however, a successful model should be able to:
  • Compute transient temperature profiles along the coal seam during gasification.
  • Determine the rate of gas and coal consumption.
  • Determine the effect of coal shrinkage or swelling on UCG performance.
  • Predict the pressure and velocity of the produced gas in a coal seam of known porosity and permeability.
  • Predict the progressive configuration of the combustion front during gasification.
  • Predict the cavity formation and dimension with time.
  • Predict the effect of operating pressure, feed inlet temperature, feed rate, feed mixture ratio, and well spacing and configuration for the gas production rate, composition, and its heating value as well as cavity dimension.
  • Predict the effect of coal bed thickness, ash content, moisture content for the gas production rate, composition, and its heating value as well as cavity dimension.
In this article, a brief overview of the UCG reactor sub-model of forward gasification is particularly discussed due to its direct impacts on coal resource recovery, cavity growth, energy efficiency and, consequently, economic feasibility.
Over the years, several approaches have been developed for the modeling of the UCG process; however, the main approaches are:
-
Packed bed model.
-
Channel model.
-
Coal slab model.
Since there are a number of mathematical/numerical models in the literature, it is most instructive to indicate the significant features and limitations of these models. Most of the earlier models were one-dimensional (1D); however, with the advancement of computational power, two-dimensional (2D) and even a few three-dimensional (3D) models were developed. In the following section, the various modeling approaches are discussed.

2.1. Packed Bed Models

The earliest models of UCG in literature include models that describe the process as a packed bed reactor. The consideration of the “in-situ” coal seam as a packed bed primarily originated from the concept of Higgins [41] considering the creation of a permeable zone between two boreholes either by reverse combustion linking (RCL) or by fracturing the coal seam using pressurized air or chemical explosives [42,43]. The resultant coal seam resembles a packed bed where coal particles are filled in the reactor. This concept is similar to the major Soviet approach to seam preparation where they included extensive drying of the seam and reverse combustion to obtain a region of enhanced permeability between the boreholes [44]. On the other hand, according to Gunn and Whitman [45], the coal seam consists of lignite or subbituminous coal can be gasified without establishing any linkage path due to their high permeability, also known as the percolation or filtration method which more closely resembles the operation of a packed bed chemical reactor. The packed bed model assumes that coal gasification occurs in highly permeable porous media with a stationary coal bed which is consumed over time [46]. Table 4 and Table 5 represent a glimpse of a list of packed bed models discussed here with their essential features and reaction rate control mechanisms of the main reactions in the gasification process, respectively.
Initially, around 1976, several 1D transient packed bed models [42,47,48] with the finite-difference approach were developed in parallel, supported by the Lawrence Livermore National Laboratory (LLNL) program. However, there is no basic difference between the models developed by Thorsness et al. [48] and Thorsness and Rozsa [42]. Thorsness and Rozsa [42] provided a detailed description of the lab-scale gasification experiment whereas Thorsness et al. [48] provided a detailed description of the development of the physical and chemical model where they made many simplifying physical and chemical assumptions, based on experimental and theoretical comparisons and correlations, which include:
-
the gas permeability.
-
the effective thermal conductivity.
-
the interphase heat-transfer coefficient.
-
the chemical reaction rates.
-
the various thermodynamic properties of each species, and.
-
the stoichiometric coefficients.
In their model, they neglected tar condensation and plugging, gas losses to surroundings, water intrusion from surroundings, heat losses, and coal bed movement due to shrinking or swelling during drying and pyrolysis in order to avoid complexity.
Winslow [47] also followed their work, except using a different numerical solution scheme. Thorsness et al. [48] considered all the main reactions for gasification (reaction 1 to reaction 10); however, for homogeneous reactions, only the water-gas shift reaction was considered by Winslow [47]. Although, for convenience, most of the researchers represent char as a molecule containing one carbon, they considered the chemical formula of char and dry coal as CHa1Ob1 and CHa2Ob2, respectively, where the composition parameters a1, b1, a2, b2 depend on the type of coal (usually obtained from the ultimate analysis). Beside the moisture content for drying, they also considered mobile condensed water which undergoes evaporation and condensation upon the removal of heating. However, for the simplicity of the calculation, they assumed a fixed velocity of mobile water. Their model considered fitted kinetic models for pyrolysis data developed by Campbell [49]. The controlling mechanisms of reaction rates that they considered are listed in Table 5. Eight gas species (N2, O2, H2O, H2, CH4, CO, CO2, tar) and five solid/condensed species (coal, char, mobile liquid water, fixed liquid water, and ash) were considered, upon which the balances were written. The reaction rates of the heterogeneous reactions (char-gas reactions) were considered to be influenced by three mechanisms: mass-transport limitation from the bulk gas to the solid surface, mass-transport limitations in any internal particle porosity, and kinetic limitations at the solid surface. This reaction mechanism was also considered for the water-gas shift reaction due to the influence of a catalyst that might be present in the coal bed. Because of the high reactivity of char as compared to pure carbon, the mass transfer limitation was considered important. However, in their model, one effective chemical rate expression (Rc) was assumed in order to account for the internal transport and surface kinetics. As a result, the total rate (RT) was considered as a kinetic internal-diffusion rate and expressed as follows:
1 R T =   1 R c +   1 R m
where Rc is linearly dependant on the mole fraction of the limiting reactant (yl) in the bulk phase and the parameter Rm is the mass transfer rate of a limiting reaction defined by:
R m = k y y l
where ky is the mass transfer coefficient.
In their models, the fluid flow was described by the Darcy equation, in which permeability was set to change with the extent of the reactions. From the experimental data of drying and pyrolysis of subbituminous Wyodak coal, they formulated a functional relationship relating permeability to porosity that was used later by a number of researchers [12,47,50] who developed their model using subbituminous coal. The relationship is as follows:
α = α 0  exp [ σ ( φ φ 0 ) ]   for φ > φ lim
where αo is the initial permeability, φo is initial porosity of coal, and φlim is the upper limit of porosity above which the permeability was assumed to be constant. Based on the experiment, the calculated value of the parameter σ was approximately 12. In their report, they considered φlim = 0.25 so that the particle size decreases with increasing porosity (φ < 0.25). As a result, the permeability was assumed to reach a constant value for porosity larger than 0.25.
Thorsness et al. [48] incorporated a quasi-steady approximation for changes associated with the gas phase that allowed them to use ordinary differential equations for gas phase mass and energy balances. The governing equations were as follows:
Mass balance:
Gas:  d u C i d x = Σ j a i j R j   ( i = 1 , 2 , .   G )
where G is the number of gas species.
It is noted that gas phase diffusion, which had been considered an important factor in the channel and coal slab model development, was neglected.
Solid:  ρ k t = M k Σ j a k j R j   ( k = 1 , 2 , 3 . K )
where K is the number of solid/condensed species.
Mobile water:  ρ w t = u w ρ w x + M w Σ j a k j R j   ( k = 1 ,   i . e . ,   o n l y   w a t e r )
Energy balance:
Gas:  d T d x =   1 u C p g [ h T ( T T s ] +   Q g
Solid:  C p s T s t =   x [ ( 1 φ ) D T s x ] +   h T ( T T s ) + Q s u w ρ w C w . T s x
Momentum balance:
The momentum balance gas phase was replaced by Darcy’s law:
p x = u α μ g
The use of Darcy’s law and some theoretical correlations they employed in their model indicate an implicit assumption of the laminar flow process.
On the other hand, Winslow [47] did not approximate any stationary phase for the basic conservation equation, although he followed the physical and chemical model of Thorsness et al. [48]. The following governing equations for mass and energy balance for both the gas and solid phases were used by Winslow [47]:
Overall gas conservation:  t ( φ ρ g ) = x ( ρ g u ) +   S m
Conservation of gas species:  t ( φ C i ) + x ( C i u ) = Σ j a i j R j ( i = 1 , 2 , .   G )
where G is the number of the gas species. Gas diffusion was neglected.
Conservation of condensed species:  t ( 1 φ ) ρ k = Σ j a k j R j ( k = 1 , 2 , 3 . K )
where K is the number of solid/condensed species.
Energy conservation of gas:  t ( φ ρ g e ) + x ( h ρ g u ) = h T ( T T s ) + Q g
Energy conservation of solid:  t [ ( 1 φ ) ρ s e s ] x [ ( 1 φ ) D T s x ] =   h T ( T T s ) + Q s
Using Darcy’s law, the gas bulk velocity is given by:
u = α μ g p x
Table 4. Some essential features of reported packed bed models.
Table 4. Some essential features of reported packed bed models.
ResearchersDimension and Time DependenceHeat TransferMass DiffusionFluid FlowThermo Mechanical FailureWater InfluxHeat Loss
ConductionConvectionRadiation
Winslow [47]1D & T D
Thorsness et al. [48]
Thorsness and Rozsa [42]1D & T D
Khadse et al. [46]1D & PS D
Uppal et al. [43]1D & PS D
Thorsness and Kang [51]1D & T D
Gunn and Whitman [45]1D & PS
Abdel-Hadi and Hsu [8]1D & T
1D = One-dimensional, D = Darcy flow, PS = Pseudo-steady-state, T = Transient.
Table 5. Reaction rate control mechanisms of the main reactions in the gasification process.
Table 5. Reaction rate control mechanisms of the main reactions in the gasification process.
ResearchersDryingPyrolysisChar ReactionWater-Gas Shift ReactionGas Phase Reaction
C + O2 → CO2C + CO2 → 2COC + H2O → CO + H2C + 2H2 → CH4CO + H2O ↔ CO2+ H2CO + 0.5O2 → CO2H2 + 0.5O2 → H2OCH4 + 2O2 → CO2 + 2H2O
Winslow [47]DPKKKKK
Thorsness et al. [48]DPKKKKKIII
Thorsness and Rozsa [42]DPKKKKKIII
Khadse et al. [46] PKKKKK
Uppal et al. [43] PKKKKKIII
Thorsness and Kang [51] KKKKKIII
Gunn and Whitman [45] ECP P
Abdel-Hadi and Hsu [8]DPPPPPE
D = Diffusion-limited, I = Infinite rate, K = Kinetic (power law) and bulk diffusion, P = Power law kinetics, EC = Experimental correlations.
The primary objective of their preliminary works was to predict and correlate reaction and thermal-front propagation rates and product-gas composition as a function of coal bed properties and process operating conditions, as well as the secondary objective of testing the correctness of the physical and chemical model they developed. Despite the simplifying assumptions as shown in Table 6, their results were generally in accordance with steady-state laboratory measurements (Table 7) Experimental run 5 was carried out using back heating in order to eliminate radial heat losses. As a result, model calculations with the assumption of no heat loss were closer to the data of experimental run 5. Although the lab-scale experiment did not represent the field-scale system as the porosity (50%) and the permeability (150 D) were too high and the particle diameter (1 cm) was too small, the agreement between the model calculations and the lab experiments indicates that the simple physical and chemical models of Thorsness et al. [48] are useful to some extent in the design of the gasification process.
Table 6. Input parameters used for model calculation and experimental runs from a 1.6 m reactor.
Table 6. Input parameters used for model calculation and experimental runs from a 1.6 m reactor.
Input ParametersCalculated (Thorsness and Rozsa [42], Thorsness et al. [48])Calculated
(Winslow [47])
Reactor experiments (Thorsness and Rozsa [42], Thorsness et al. [48])
Run 4 (no backheating)Run 5 (with backheating)
Coal DimensionL = 150 cm,
D =15 cm
L = 150 cm,
D =15 cm
L = 150 cm,
D =15 cm
L = 150 cm,
D = 15 cm
Initial Porosity50%44.4%50%50%
Initial Permeability150 Darcy 150 Darcy150 Darcy
Coal Particle Size1 cm1 cm1 cm1 cm
Injected Gas Flow Rate2 × 10−4 gmole/cm2·s2.1 × 10−4 gmole/cm2·s2 × 10−4 gmole/cm2·s2 × 10−4 gmole/cm2·s
Initial Feed Temperature450 K450 K450 K450 K
Pressure4.8 Atm4.8 Atm4.8 Atm4.8 Atm
Injected Gas Composition
O214.1%13.4%13.4%13.4%
N224.0%22.0%22.0%22.0%
H2O83.5%84.4%84.4%84.4%
Steam/O2 Ratio6666
Table 7. Steady-state temperature and major constituents of outlet gas composition.
Table 7. Steady-state temperature and major constituents of outlet gas composition.
ParametersCalculated (Thorsness and Rozsa [42], Thorsness et al. [48])Calculated (Winslow [47])Measured (average over run) (Thorsness and Rozsa [42], Thorsness et al. [48])
Run 4 (no back heating)Run 5 (with back heating)
Peak Temperature (K) 1250 1200
Speed of burn front18 cm/h18 cm/h 18 cm/h
Product Gas Composition (dry basis)
H243.2%43.8%44.9%44.6%
CH46.6%4.0%6.4%6.9%
CO22.0%22.0%16.5%19.0%
CO228.2%30.2%33.2%29.5%
Beside the model validation with experimental data, they also illustrated various features of their simulations with time or spatial variations and compared them with the experimental results where possible. The solid and gas temperature distribution along the length in their system at different times indicated the propagation of the reaction and drying/pyrolysis fronts at different speeds (see Figure 3).
Figure 3. Calculated gas and solid temperature at three different times for the 1.5 m reactor simulation [48].
Figure 3. Calculated gas and solid temperature at three different times for the 1.5 m reactor simulation [48].
Energies 08 12331 g003
The different fronts that were growing at different speeds associated with the movement of the gasification process were also recognized from the plot of the reaction rate as a function of bed length (Figure 4). Those were the coal drying front, then the coal pyrolysis front, followed by the char reaction front considering moving upstream from the outlet. The double-peak of the water-gas shift reaction was due to the changes in equilibrium limits because of changes in temperature (Figure 4). The gas concentration along the bed length showed that the inlet gas composition remained unreacted until the reaction front (Figure 5). However, after an initial rapid change followed by a more gradual one, a sharp change in gas concentration was observed due to the drying-pyrolysis front. Winslow [47] also observed the same trend.
Figure 4. Distribution of principal reaction rates at t = 1 h 13 min [47].
Figure 4. Distribution of principal reaction rates at t = 1 h 13 min [47].
Energies 08 12331 g004
Figure 5. Calculated gas phase concentration at time = 1.4 × 104 s for the 1.5 m reactor simulation [48].
Figure 5. Calculated gas phase concentration at time = 1.4 × 104 s for the 1.5 m reactor simulation [48].
Energies 08 12331 g005
A comparison of the position of the calculated reaction and drying front with experimental data showed that the results were in good agreement with the reaction front; however, the predicted drying front movement rate was found higher than the experimental findings (see Figure 6). According to them, the latter observation was due to some uncertainty of the product gas mix and the heat losses occurring during the experiment at the drying zone.
Figure 6. Comparison of calculated and experimental reaction and drying front movements [48].
Figure 6. Comparison of calculated and experimental reaction and drying front movements [48].
Energies 08 12331 g006
Thorsness et al. [48] extended their works by reporting very limited parametric study. For example, with increasing oxygen concentration in the feed gas (holding a constant gas flow rate), they observed an increased gas/solid temperature and the drying rate of the reaction front moved linearly with oxygen concentration. An increased ratio of CO to CO2 resulting from the increased temperature generation in the system was also detected. With the increase of total dry gas production, a consistent CO/CO2 ratio was also noticed with the higher oxygen feed. However, for energy recovery per mole of O2, they indicated an optimum value of oxygen consumption at a mole fraction of approximately 0.2 (see Figure 7). This observation is actually set a limit for the steam-to-oxygen ratio at a value of 4 for maximum energy recovery.
Figure 7. Calculated energy recovery for three levels of oxygen feed concentration [48].
Figure 7. Calculated energy recovery for three levels of oxygen feed concentration [48].
Energies 08 12331 g007
There was an attempt by Thorsness et al. [48] to simulate a field-scale system with length equal to 15 m, porosity equal to 0.05, permeability equal to 20 D, outlet pressure equal to 1 atm, and particle diameter equal to 20 cm. However, they could not validate the changes of the process variable except the developing temperature profile that was similar to their lab-scale simulation.
Although they showed time and distance changes in the process variable in the principal direction of gas flow, they did not conduct a detailed sensitivity analysis to investigate the effect of coal particle diameter, porosity, permeability, injection gas flow rate, and other parameters in their model despite the fact that the basis of the model formulation was a 1 cm particle diameter. This imposes the use of their model for specific conditions. The applicability of their model to field conditions where the particle diameter is approximately equal to 20 cm is limited because of the lack of sensitivity of the results to particle size and the validation of field trials. The particle size is particularly important because, in their physical and chemical model, the char reaction rates and other transport coefficients such as heat transfer and mass transfer coefficients are directly or indirectly related to the particle sizes. In addition, because of high initial porosity (φ > 0.25), they considered a constant permeability in their simulation. There apparently was no use of their formulated functional relationship relating permeability to porosity for considering high initial porosity. The shape of the combustion zone and, hence, the actual gasification volume of the coal could not be determined since the model was only one-dimensional.
Because of its transient nature, the model developed by Winslow [47] can be applied for the cases where rapidly changing transient behaviour is important. However, no stark difference between the model outcomes of transient and quasi-steady-state models was observed, and that is probably why some researchers adopted the quasi-steady-state consideration for ease of calculations.
After three decades, in 2006, Khadse et al. [46] developed a simple 1D model following the physical and chemical model as well as the pseudo-steady-state fluid flow model described by Thorsness et al. [48]. However, their model differs from the model developed by Thorsness et al. [48] in that the drying and mobile water evaporation/condensation reactions were not considered and only coal and char were considered as solid phase. Nevertheless, their work gets attention as they analyzed the effects of various operating and model parameters on the temperature and gas phase and solid compositions in UCG that were not completely incorporated in the models discussed above. Their input parameters for the base case were slightly different than that of Thorsness et al. [48] as can be seen from Table 8. However, at least for one run, they simulated the experimental conditions of Thorsness et al. [43] in order to compare their calculated exhaust gas composition with the experimental data (Figure 8). However, they could not compare the calculated results of Thorsness et al. [48] that were very close to the experimental results quantitatively, although they followed the development of the model of Thorsness et al. [48]. This is possibly due to considering certain parameters that lack specific information in the literature model and neglecting the drying and mobile water evaporation/condensation reactions that Thorsness et al. [48] considered. However, their result can be considered only in a qualitative agreement with the experiments of Thorsness et al. [48].
Table 8. Comparison of input parameters of Khadse et al. [46] with experimental run 5 of Thorsness et al. [48].
Table 8. Comparison of input parameters of Khadse et al. [46] with experimental run 5 of Thorsness et al. [48].
ParametersKhadse et al. [46] (Base Case)Experimental Run 5 [48]
Input Parameters
Coal DimensionL = 100 cm, D = 15 cmL = 150 cm, D = 15 cm
Initial Porosity20%50%
Coal Particle Size1 cm1 cm
Injected Gas Flow Rate2 × 10−4 gmole/cm2·s2.1 × 10−4 gmole/cm2·s
Initial Feed Temperature430 K450 K
Pressure4.8 Atm4.8 Atm
Injection Gas Composition
O215.4%13.4%
N27.6%22.0%
H2O77.0%84.4%
Steam/O2 Ratio56
Figure 8. Comparison of calculated exhaust gas compositions with experimental results from Thorsness et al. [48]; (a) Khadse et al. [46]; (b) Thorsness et al. [48]).
Figure 8. Comparison of calculated exhaust gas compositions with experimental results from Thorsness et al. [48]; (a) Khadse et al. [46]; (b) Thorsness et al. [48]).
Energies 08 12331 g008
The simplicity of the model enabled authors to investigate the effects of various parameters such as O2 concentration, steam/O2 ratio, and inlet pressure in their model. They found that increasing oxygen concentration (holding the steam fraction constant) in the feed increases the propagation rate of the reaction front and the temperature in the reaction zone. Similarly, increasing the steam fraction (holding the oxygen concentration constant) increased the propagation rate of the reaction front; however, a decrease in temperature was observed. In the absence of nitrogen, they found an optimum outlet product with a steam/oxygen ratio of 1.5. They also found that input parameters did not influence the pyrolysis reaction. They also reported the variation of the gross calorific value of the exhaust gas, on a dry basis, at different pressures with time (Figure 9) where the gross calorific value was calculated using the following formula obtained from Harker and Backhurst [52]:
C a l o r i f i c   V a l u e   ( k J m o l e ) = H i × y i
where Hi is the heat of the combustion of the gas “i” (kJ/mol), and yi is the mole fractions of the gas “i” in the product gas on a dry basis. However, the calculation was based on the exit concentration of CO, H2, and CH4 only. Pressure was found to have no effect on the process (Figure 9). In all the cases, a steady flow rate and composition in the outlet was reached after 5000 s of process simulation.
Figure 9. Gas calorific value at different total pressures with time [46].
Figure 9. Gas calorific value at different total pressures with time [46].
Energies 08 12331 g009
Recently, Uppal et al. [43] also developed another 1D packed bed model adopting the existing model of Thorsness et al. [48] with some modifications in the model structure and solution strategy. In order to observe their model capability, they considered a bed length of 4.8 m only and observed various features of their simulations with qualitative agreement with the literature. However, for validating their model with the experimental data, they considered the experimental data from their pilot plant which had dimensions of 150 m × 125 m that contained an array (7 × 6) of wells connected by a network of pipes. For comparing with the model calculation, only two consecutive wells that were 25 m apart from each other, were considered. Table 9 shows their input parameters, some of which are significantly different from the model that they followed.
Table 9. Input parameters for UCG packed bed model developed by Uppal et al. [43].
Table 9. Input parameters for UCG packed bed model developed by Uppal et al. [43].
ParametersValue
Length of Reactor25 m
Coal TypeLignite B
Injected GasAir
Injected Gas Flow Rate2 × 10−4 gmole/cm2·s (for observation of model capabilities)
Initial Feed Temperature430 K
Pressure6 Atm
Initial Coal Density1.25 g/cm3
In their pilot plant, the injected gas (air only) flow rate was varied with time (Figure 10) and the resultant heating value and the exhaust gas composition were recorded. In their model, they calculated the exit gas heating value and composition utilizing the experimental inlet gas flow rate. In addition, they used a nonlinear optimization technique in order to compensate the uncertainty in coal and char by optimizing the composition parameters of some pyrolyzed products (H2, CH4, and char). The calculated results were in a reasonable agreement with the experimental results (Figure 11a,b). However, according to them, a better prediction can be obtained by increasing the optimization variables. Their model lacks some detailed information such as particle diameter, porosity, permeability, etc., which are considered important for understanding a UCG process. Nevertheless, their model indicates that a control strategy can be employed by manipulating the inlet flow rate to maintain a desired heating value in the presence of disturbances.
Figure 10. Flow rate of the injected air with time [43].
Figure 10. Flow rate of the injected air with time [43].
Energies 08 12331 g010
Figure 11. Experimental and simulated (a) heating values and (b) mole fractions of important gases of the product gas mixture with time (reproduced from Uppal et al. [43]).
Figure 11. Experimental and simulated (a) heating values and (b) mole fractions of important gases of the product gas mixture with time (reproduced from Uppal et al. [43]).
Energies 08 12331 g011
In 1986, Thorsness and Kang [51] developed a generalized 2D model for describing reacting flows through the packed bed using the following governing equations:
Mass balance:
Overall gas conversion:
( φ C ) t =     . ( v ¯ φ C ) + Σ j a i j R j + S m ( i = 1 , 2 , .   G )
Conservation of gas species:
t   ( φ C i ) =   . ( v ¯ . φ C i ) + . [ C D i y i ] + Σ j a i j R j +   S m i
Overall solid conservation:
t   [ ρ s ( 1 φ ) ] =   . [ ( 1 φ ) ρ s v s ¯ ] +   Σ S s k
Solid species conservation:
t   [ ρ s Y s k ( 1 φ ) ] =   . [ ρ s v s ¯ Y s k ( 1 φ e ) ] +   S s
Energy balance:
They assumed an identical gas and solid temperature and incorporated one energy equation for the combination of gas and solid as follows:
t   [ φ Σ ( C i h i ) + ( 1 φ ) ρ s Σ ( Y s k h k ) ] =   . [ Σ h i J l ¯ ]   . [ ρ s v s ¯ Σ ( Y s k h k ) ] + . [ k e f f T ] + Σ r a t e i H i
Considering the negligible difference between the gas and solid temperature obtained by the earlier models [47,48], the assumption of identical gas and solid temperature seems to be justified.
In their model, they incorporated diffusion effects, wall transport, and also char combustion and water-gas shift reaction rates based on Shell Progressive (SP) and Ash Segregation (AS) reaction models. In the SP model, a core of unreacted solid was assumed to be surrounded by a shell of ash through which the gas phase reactants diffuse (Figure 12a). On the other hand, the AS model assumed continuous exposure of unreacted material to the gas stream due to the ash segregation (Figure 12b).
Figure 12. Different reaction regimes in packed bed model.
Figure 12. Different reaction regimes in packed bed model.
Energies 08 12331 g012
Reaction rates for these two models were as follows:
1 R S P =   1 R b u l k   m a s s   t r a n s f e r +   1 R a s h - l a y e r   d i f f u s i o n +   1 R d i f f u s i o n   i n t o   +   n t r i n c i s   r e a c t i o n 1 R A S =   1 R b u l k   m a s s   t r a n s f e r +   1 R d i f f u s i o n   i n t o   +   n t r i n c i s   r e a c t i o n
Because of the generalized nature of their model, they tested various cases (i.e., transient concentration and thermal wave, steady heat transfer phenomena, etc.) for non-reacting packing materials through the packed bed and compared them favourably with other available studies. Although their generalized model is 2D, only one case (steady heat transfer phenomena) was solved using the 2D model. For all other cases, the 1D model was considered. For validation of a UCG model, they calculated gas composition, temperature, and carbon fraction considering a steady 1D model with very limited gas species and compared the results with the analytical data obtained from literature. They analyzed the following three disparate situations related to the UCG process which were not considered in the earlier models:
-
Case 1: Transient water injection into the midpoint of the packed bed during gasification.
-
Case 2: Coal wall drying using a hot gas flow.
-
Case 3: Wall regression during gasification.
Table 10 shows the input parameters for the above situations. For the first case, seven gas species (N2, O2, H2, CO, CO2, water vapor, and CH4) and eight reactions (reactions 3–11) were considered. However, drying and pyrolysis reactions were ignored. During gasification process, when the product gases reached a reasonably steady state, water was injected at the midpoint of the bed, slightly downstream of the reaction front, and continued until the mole fraction of CO and CO2 were observed to approach steady values. However, the liquid water influx was assumed to turn into steam instantly with a flow rate of 8 × 10−4 mol/s. A sudden change in CO and a gradual change in H2 were noted after the water injection; however, the gases undergoing changes returned to pre-injected conditions as seen in Figure 13.
Table 10. Input parameters for the three situations considered by Thorsness and Kang [51].
Table 10. Input parameters for the three situations considered by Thorsness and Kang [51].
Input ParametersCase 1Case 2Case 3
Packed Bed DimensionL = 150 cm, D = 5 cmL = 25 cm, D = 10 cmL = 1 m, D = 1 m
Initial Gas Temperature900 K900 K900 K
Initial Bed Temperature 300
Gas Composition
H2O74%100% (or N2)66.6%
O226% 33.3%
Injected Gas Flow Rate1.5 × 10−3 mole/s8 × 10−3 mole/s
(1 mol/s-m2 of bed)
6 moles/s
Figure 13. Product gas changes during gasification with midpoint water injection [51].
Figure 13. Product gas changes during gasification with midpoint water injection [51].
Energies 08 12331 g013
For case 2, a uniformly permeable non-reacting bed with a water-saturated wall at steam temperature was assumed through which hot gas was passed. It was observed that the time to reach steady state increased rapidly with a decreasing flow rate as if the average wall-drying rate is a function of the injected gas flow (Figure 14). However, the wall-drying rate was found to be limited by the available energy for too low injected hot gas due to the insufficiency of the total heat injected.
Figure 14. Drying rate vs. flow rate (log-log plot, from Thorsness and Kang [51]).
Figure 14. Drying rate vs. flow rate (log-log plot, from Thorsness and Kang [51]).
Energies 08 12331 g014
For case 3, the model system was considered to be filled with rubble material consisting of ash in the center and char near the walls and at the top. The walls were considered as coal that can proceed with gasification reactions. The objective of this model system was to determine the thickness of the char bed at the wall that would lead to a self-sustaining system. Calculations were performed using all seven gas species. The computed wall regression rate (i.e., the carbon production rate) was found nearly constant with the increasing char layer thickness, while the carbon consumption varied linearly with the char thickness (Figure 15). A wall char bed thickness of 7.5 cm was found close to a self-sustaining system, since the rate of carbon consumption was nearly equal to the carbon production. However, the calculated average wall regression rate of 7.7 × 10−7 m/s (0.07 m/day), when comparing the characteristic field result (~0.5 m/day), indicates that their model is missing some additional physics.
Figure 15. Carbon production rate vs. wall char layer thickness [51].
Figure 15. Carbon production rate vs. wall char layer thickness [51].
Energies 08 12331 g015
In 1976, in parallel to LLNL, Gunn and Whitman [45] also developed a 1D linear mathematical model separately for UCG. Their work is particularly important as their model is the only packed bed model cited in the literature that was compared with the field trials (Hanna, Wyoming) and exhibited a fairly good agreement with the experimental data (Table 11). However, field test data was considered at a time span when a steady-state condition was observed in order to conform closely to the assumption of their model. In addition to the gas composition, temperature, and heating value, they also determined the ratio of the gas production rate to the gas injection rate, the effect of reservoir water influx on the heating value of gas, the distribution of thermal energy and thermal efficiency during the gasification process that closely matched the experimental observations and, thus, provided enhanced understanding of the UCG process.
Table 11. Comparison of experimental and calculated data for an air injection rate = 1631 Mcf/day [45].
Table 11. Comparison of experimental and calculated data for an air injection rate = 1631 Mcf/day [45].
ComponentsComposition (Mole Percent)
Field Test ResultsCalculated Data (Average Value)
H218.5718.60
CH44.104.92
N2 + Ar47.8147.13
CO16.3515.83
CO212.3313.29
H2S0.040.07
Ethan +0.800.14
Gas Production Rate (Mcf/day)27282734
Gas Production Temperature, °F493 (measured at the surface)533 (measured at the bottom of the well)
Maximum Temperature, °F18711951
Heating Value, Btu/scf170.6164.6
Thermal Efficiency (%)89.287.4
Unlike the above models, the water-gas shift reaction and methanation were not included in their model. Methane was considered to be produced solely through pyrolysis. Only two reactions, water-gas and combustion reactions, were considered. The reaction rate of the water-gas reaction was adopted from the work of Gadsby et al. [53]; however, the constants were determined from the gasification of the pitch coke. For the combustion reaction, the reaction data of Lewis et al. [54] was used that was obtained from the combustion of metallurgical coke. Even the correlation of the thermal conductivity they used was not intended for subbituminous coal. Drying was not considered in their model. Despite this and the large errors in the reaction rate, they observed only negligible errors in the model predictions.
In their model, water vapor, carbon monoxide, hydrogen, carbon dioxide, oxygen, inert gases (nitrogen and argon), and devolatilized material from the coal were considered. For convenience, the volatile products resulting from the pyrolysis of coal were treated as a single pseudocomponent; however, the pseudocomponent was broken down into individual gases for the purpose of calculating total gas composition. From the proximate and ultimate analysis of the experimental subbituminous coal (from the Hanna field test), they also formulated a correlation between the weight fraction of volatile matters and the temperature with high accuracy. In order to obtain the composition of devolatilized gas, they considered the experimental devolatilized data of the Hanna I coal seam at 900 °C along with the assumption that the tars (16.2%) and light oil (0.9%) and water (13.3%) in that analyses cracked to the stoichiometric proportion of methane, carbon monoxide, and hydrogen. Their assumption led to the composition for devolatilized products listed in Table 12.
Table 12. Composition of devolatilized products assumed by Gunn and Whitman [45].
Table 12. Composition of devolatilized products assumed by Gunn and Whitman [45].
ComponentsPercent
H242.26
CO28.53
CO23.79
CH424.36
H2S0.37
Ethane +0.69
Total100
Average molecular wt.14.76
Their model transformed the set of partial differential equations into a set of ordinary differential equations by considering the pseudo-steady-state approximation. They considered nine integrated continuity equations (two for overall solid and mass balances and others for individual gas species), a partially integrated energy equation, and a macroscopic material balance equation from which a constant combustion front velocity could be calculated iteratively. The assumption of the constant combustion velocity seems to be valid, according to Figure 6, with what was observed by Thorsness et al. [48].
Their model predicted the exhaust dry gas composition and obtained a close agreement with the experimental gas composition history with the exception of methane (Figure 16). However, the deviation of methane was believed to be a result of considering a high concentration of methane in the devolatilized products. The only parameter fitted in their model from the field trial was the water influx in order to get a good prediction for the data that interferes with the water influx rate. Their calculated results were also closely matched with experimental data for ratio of gas production/air injection and gross heating value (Figure 17 and Figure 18). However, the declination of those quantities with time was believed to be the result of the gradual increase of water flux during field test. They calculated the heating value of gas by changing the ratio of water influx to air injection rate and observed a maximum heating value at a ratio of about 0.15. Usually the maximum heating value at a particular ratio of water intrusion to air injection rate implies the control strategy of closely maintaining the optimum value either by controlling the water intrusion rate by pressure of any other means or by changing the air injection rate.
Figure 16. Comparison of calculated exhaust gas compositions from Gunn and Whitman [45] with experimental results from the Hanna field trial.
Figure 16. Comparison of calculated exhaust gas compositions from Gunn and Whitman [45] with experimental results from the Hanna field trial.
Energies 08 12331 g016
Figure 17. Comparison of calculated data and field test: (a) gas production/air injection ratio and (b) gas heating value [45].
Figure 17. Comparison of calculated data and field test: (a) gas production/air injection ratio and (b) gas heating value [45].
Energies 08 12331 g017
Figure 18. Effect of reservoir water influx on heating values of gas [45].
Figure 18. Effect of reservoir water influx on heating values of gas [45].
Energies 08 12331 g018
The calculated values of the distribution of thermal energy released by coal combustion were also in good agreement with the experimental values (Table 13). From the experimental data, it was observed that only a small amount of energy from a 30-foot coal seam was lost to the overburden and base rock. This implies that neglecting heat losses to the surrounding rocks does not hamper a model outcome too much. The last column of Table 13 shows a higher heating value of the product gas when injected air was preheated at 480 °F. In order to get higher economic value and energy conservation, they suggested using a part of the heat of vaporization of the produced water to heat the inlet gas. They also reported that in all cases, the reaction zone was very narrow, generally two feet or less. This fact was verified by both model calculations and by temperature measurements in the observation wells.
Table 13. Distribution of thermal energy from in situ coal combustion [45].
Table 13. Distribution of thermal energy from in situ coal combustion [45].
Calculated, %Field Test, %Calculated with Preheating, %
Produced Gas
Gross heating value87.489.191.3
Sensible heat....5.85.22.9
Heat of vaporization of water vapor in produced gas3.43.32.4
Heat Loss to Surrounding Rock1.33.81.3
Heat Loss in Ash2.10.62.1
Total100100100
Although they obtained a good agreement with the experimental results for a number of phenomena, they could not quantitatively predict the variation of the temperature and mole fraction along the bed length with the experimental value because of the poor consideration of the thermal conductivity, heat capacities, and the uncertain accuracy of the reaction rates. However, they obtained a qualitative agreement with the experimental values. Despite very close agreement with the experimental data, their model is not sufficient for predicting all aspects of gasification as their model neglected some important phenomena, such as drying, CO2 gasification, the water-gas shift reaction, gas phase reactions, etc. However, because of the low water/air ratios that were experienced with the Hanna tests, there might be little or no water-gas shift reactions. As a result, the Gunn model was not hampered by neglecting this reaction. In their model, there is also a lack of sensitivity analysis for the kinetic parameters they used. In addition, they included water influx data from a particular field trial. However, water influx is a highly variable phenomenon that depends on coal seam properties and hydrology. All these lacking parameters limit their model for general applicability in the UCG process.
Abdel-Hadi and Hsu [8] extended previous models by developing pseudo-2D geometry with a moving burn front in the axial direction. A rectangular domain with a length of 1.5 m and width of 1 m was used in their model. Their governing equations are nearly similar to the equation considered by Winslow [47] and Thorsness and Kang [51]; however, they included carbon consumption in the reaction zone in order to track the burn front, and the equation is given by:
ρ s w c v c = 12   0 x r z ( Σ j R c i j ) d x
where wc is the carbon mass fraction, vc is the velocity of the combustion front, and xrz is the total length of the reaction zone. They also performed an immobilization transformation of coordinates in order to formulate the irregular front motion as follows:
ξ   = x X ( y , t ) H 1 X ( y , t ) η = y H 2
All other governing equations were transformed to these new coordinates (Figure 19). This allowed them to calculate the progressive configuration of the gasified zone at various stages of gasification (Figure 20) which, in turn, facilitated the calculation of the rate of coal seam consumption (Figure 21). The conversion rate of the coal seam is found fairly constant. In order to gain confidence for this model, they have compared their model with the laboratory results reported by Thorsness and Rozsa [42] and obtained a good agreement with the experimental data.
Figure 19. Schematic description of the problem: (a) conventional coordinate system; (b) immobilization coordinate system (redrawn from Abdel-Hadi and Hsu [8]).
Figure 19. Schematic description of the problem: (a) conventional coordinate system; (b) immobilization coordinate system (redrawn from Abdel-Hadi and Hsu [8]).
Energies 08 12331 g019
Figure 20. Progressive configurations of combustion fronts during gasification [8].
Figure 20. Progressive configurations of combustion fronts during gasification [8].
Energies 08 12331 g020
Figure 21. Consumption of coal during gasification [8].
Figure 21. Consumption of coal during gasification [8].
Energies 08 12331 g021
The packed bed models have been validated with laboratory experiments to some extent. These models exhibit good agreement for the gas composition. In order to calculate the heat recovery and gas composition, the model can be very effective. However, the extension of these models to field trials is infeasible, since other cavity growth mechanisms such as thermo-mechanical failure could not be incorporated into the models. Also, as pointed out by Winslow [47], this method requires a fine grid at the vicinity of the reaction front that limits its applicability to three-dimensional field-scale trials. However, recent advancement of computational power can easily overcome this limitation.

2.2. Channel Models

Beside packed bed models, a number of channel models have also been developed in the first decades of modeling. The channel model assumes that coal is gasified only at the perimeter of the expanding permeable channel [19]. In this approach, the UCG process is represented by an expanding channel when two distinct zones of rubble/char and open channel exist. This approach is considered due to the observation of the formation of the open channel structure after the gasification phase is terminated in different field tests of coal seams [55,56]. Figure 22 shows the basic concept and physics behind this approach: “Air or oxygen flows down the central channel and is convected by turbulent flow to the boundary layer along the channel wall. The oxygen diffuses through the boundary layer to the solid surface and reacts. The hot combustion gases diffuse back through the boundary layer to the channel” [40]. The channel model is more useful for analyzing sweep efficiency [19].
Table 14 and Table 15 represent a glimpse of a list of channel models discussed here with their essential features and reaction rate control mechanisms of the main reactions in the gasification process, respectively.
Figure 22. Reactions and transport phenomena in channel model (redrawn from Gunn and Krantz [40]).
Figure 22. Reactions and transport phenomena in channel model (redrawn from Gunn and Krantz [40]).
Energies 08 12331 g022
Table 14. Some essential features of reported channel models.
Table 14. Some essential features of reported channel models.
ResearchersDimension and Time DependenceHeat TransferMass DiffusionFluid FlowThermo Mechanical FailureWater InfluxHeat Loss
ConductionConvectionRadiation
Dinsmoor et al. [44]2D and T P
Eddy and Schwarz [57]1D and TM
Luo et al. [58]2D and T
Batenburg et al. [55]1D and SS D
Pirlot et al. [59]2D and S D
Kuyper [56,60]2D and TNS
Perkins and Sahajwalla [61]2D and TNS
D = Darcy flow, P = Plug flow, M = Mixed flow, NS = Navier-Stokes, S = Steady State, SS = Semi-Steady State, PS = Pseudo-Steady State, T = Transient.
Table 15. Reaction rate control mechanisms of the main reactions in the gasification process.
Table 15. Reaction rate control mechanisms of the main reactions in the gasification process.
ResearchersDryingPyrolysisChar ReactionWater-Gas Shift ReactionGas Phase Reaction
C + O2 → CO2C + CO2 → 2COC + H2O → CO+ H2C + 2H2 → CH4CO + H2O↔ CO2+ H2CO + 0.5O2 → CO2H2 + 0.5O2 → H2OCH4+2O2→ CO2+2H2O
Dinsmoor et al. [44] KKK PP
Luo et al. [58] PPPPPPPP
Batenburg et al. [55] EEEE E
Pirlot et al. [59] EEEEEEE
Kuyper [56,60] D M
Perkins and Sahajwalla [61] P M
D = Diffusion limited, E = Equilibrium, K = Kinetic (power law) and bulk diffusion, M = Turbulent mixing limited, P = Power law kinetics.
Dinsmoor et al. [44] developed a steady-state, 1D channel model by assuming that the gasifier behaves as an expanding cylindrical cavity in the coal seam with reactions taking place at the walls. In their model, for simplicity, no pyrolysis reactions were considered. Heat transfer included conduction for solids only; however, both convection and radiation were included between the wall and gas. Axial heat conduction in the gas phase was neglected. They also considered water influx as evenly distributed along the length of the tube. Char reactions (reactions 3–5) and two gas phase reactions (reactions 9 and 10) were considered in their model. Because of slow channel evolution, they incorporated a pseudo-steady-state assumption for changes associated with the gas phase. They treated a fully developed flow process with the channel despite the initial cylindrical channel diameter (0.3 m) being much smaller than the channel length (60 m). Considering forced convection as the dominating mechanism of mass transfer, they simulated coal gasification with a forced convection mass transfer correlation. For heterogeneous reaction kinetics, they considered the surface reaction rate constant and wall diffusion resistances. For wall diffusion resistances, the mass transfer coefficients were calculated from a standard correlation for the turbulent flow of gases in tubes. Inlet gas was assumed to be air at 330 K and 6.9 atm with no water added. Their predicted gas compositions and temperature profile showed clear evidence of the presence of a separate oxidation and reduction zone; however, a very high temperature (2400 K) was noted near the inlet (Figure 23). According to them, this high temperature was due to neglecting heat loss at the end of the tube and including the radiation heat transfer. However, it appears that apart from their explanations, neglecting the pyrolysis reaction could be another reason for getting high temperature initially, as during ignition, a significant amount of heat is used to create char by the pyrolysis of coal.
Figure 23. Typical channel temperature (gas and wall) and gas composition profiles without water influx [44].
Figure 23. Typical channel temperature (gas and wall) and gas composition profiles without water influx [44].
Energies 08 12331 g023
The quality of the predicted product gas observed was inferior to the gas quality usually obtained in the packed bed model for similar situations. For a constant blast velocity, the reaction rates were observed to decrease with the evolution of the channel due to the constant mole fraction of the oxygen. As a result, the oxidation zone became longer which, in turn, was responsible for increased heat losses and the deterioration of gas quality. Compared to the packed bed combustion front (approximately 0.2 m), the oxidation zone (approximately 20 m) was much longer. However, these observations were not supported by any field observations. Apparently, because of the long oxidation and reduction zones, they concluded that “a successful UCG system cannot be operated in the channel regime”.
Almost at the same time, the above conclusion was negated by the work of Schwartz et al. [62] as they found an increase of mass transfer by several orders of magnitude when natural convection in channels was considered as the controlling mechanism of mass transfer instead of forced convection alone. The intense mixing of the injected blast with the boundary layer gases depends on the Grashof number (Gr) of the natural convection flow. A good mixing begins when the Grashof number exceeds 108 and becomes intense at about 1010 [63]. However, for dominant natural convection flow, the gases within the boundary layer produce a circulatory flow within the cylinder and a good mixing is obtained. However, Schwartz et al. [62] assumed two different modes of the circulatory pattern based on coal seam thickness as depicted in Figure 24. Figure 25 shows the difference of cavity growth with and without considering the natural convection. Later, it was established by mathematical correlation that the models assuming forced convection mass transfer as the controlling mechanism for cavity growth were not born out of field trials [15].
Figure 24. Assumed circulation patterns in the thick seam and thin seam configurations (redrawn from Schwartz et al. [62]).
Figure 24. Assumed circulation patterns in the thick seam and thin seam configurations (redrawn from Schwartz et al. [62]).
Energies 08 12331 g024
Figure 25. A comparison of cavity growth via forced convection and natural convection [62].
Figure 25. A comparison of cavity growth via forced convection and natural convection [62].
Energies 08 12331 g025
Schwartz et al. [62] were the first investigators to consider the natural convection as the controlling mechanism of heat and mass transfer in UCG cavities. In their later paper, Eddy and Schwarz [57] developed a 2D model and described the evolution of the cavity based on the movement of the cavity wall. The blast from the injection well was transported through the reacting walls of a constant temperature via the convection process and ultimately diffused to the wall and reacted with the coal. The oxygen diffusion rate was calculated by natural and forced convective heat and the mass coefficient. The roof and floor of the cavity were assumed to be non-reacting surfaces. They considered neither subsidence nor rubblization. They included the water-gas shift reaction and two other gas phase oxidations along with three char reactions (reactions 3–5) as considered by Schwartz et al. [62]. In order to capture the experimentally proven teardrop shape of the cavity formation, they divided the cavity growth into different regions, i.e., a hemispherical region in the vicinity of the injection well and a series of non-equal diameter cylinders (step cylinder) downstream from the injection well, according to the flow process and a cylindrical link zone (Figure 26). At any location, the cavity was assumed to be a cylindrical cross-section until the cavity diameter was equal to the seam thickness, which is when the model switches to a rectangular cross-section. It was assumed that all coal volume located in the hemispherical region was consumed in the reaction with oxygen. The volume of char was calculated from the mass of char burned due to the combustion and the density of the char. The half-body volume was then fit to the calculated char volume. Flow was treated as an entrance region of pipe, as the calculated entrance length was much larger than the vertical well distance. The entrance length was calculated using the following equation:
L e d = 0.0288 R e d
where Le is the entrance length and d is the pipe diameter. The minimum Reynolds number in the system was 4000, which resulted in an entrance length of 48.8 m. The entrance length was considered to increase with the increase of the cavity diameter. According to them, as fresh blast proceeds down the axis and products of combustion are released at the walls, the forced and natural convection results in a swirling flow that ultimately helps mix the flows and enhance the flow of oxygen to the vicinity of the wall where boundary layer convective diffusion dominates. The transport of heat and species inside the cavity was considered to be governed by empirical correlations for the turbulent transport phenomena in enclosures. The model was able to reproduce the results of the Hanna II and Pricetown field trials qualitatively.
Figure 26. Conceptual sketch of the linked vertical well computer model geometry (redrawn from Eddy and Schwarz [57].
Figure 26. Conceptual sketch of the linked vertical well computer model geometry (redrawn from Eddy and Schwarz [57].
Energies 08 12331 g026
Recently, Luo et al. [58] extended the Schwartz et al. [57] model by including heat transfer and more coal wall and gas phase reactions (reactions 3–10). Flow inside the cavity was solved based on irrotational fluid flow inside an enclosure, which describes velocity potential based on geometric features of the enclosure. The geometry is shown in Figure 27. The stream function and velocity potential were calculated using the following equations:
Steam function:
Ψ =   m ˙ 4 π cos ( π θ ) +   U r 2 2 sin 2 ( π θ )
Velocity potential:
v =   m ˙ 4 π r +   U r cos   ( π   θ )
where U is the velocity of the uniform stream and is the mass flow rate. Pressure distribution in the cavity was determined using Bernouli’s equation:
p + 1 2   ρ u 2 =   p +   1 2   ρ U 2
Figure 27. UCG cavity defined in the work of Luo et al. [58].
Figure 27. UCG cavity defined in the work of Luo et al. [58].
Energies 08 12331 g027
In the Luo model, Fluent was used to predict the cavity shape at different times based on heterogeneous reaction rates at the cavity wall. In order to quantify the turbulent intensity, they used the standard k-ε model for the turbulent kinetic energy and dissipation. For heterogeneous and homogeneous reactions, they used kinetics/the diffusion control process and a finite rate/eddy dissipation model, respectively. Their model was validated with Chinchilla field trials. The calculated coal consumption closely matched the Chinchilla trial field data (Figure 28). For sweep cavity geometry, their model was further validated with the data from the Hanna II and III UCG trials. There is less than 5% error between the results generated with the model (Hanna II: 2436.6 tons of coal gasified in 25 days; Hanna III: 4139.3 tons in 38 days) and reported data (Hanna II: 2500 tons of coal gasified in 25 days; Hanna III: 4200 tons in 38 days). Their model predicts a hemispherical shape for the cavity geometry; however, the model is limited since heat and mass transfer characteristics of the cavity are unknown. Also, the coupling of this model with the mechanical failure of the coal would be cumbersome, since the accumulated rubble in the cavity changes the transport phenomena inside the cavity.
Figure 28. Coal consumptions for modeling and Chinchilla trial data [58].
Figure 28. Coal consumptions for modeling and Chinchilla trial data [58].
Energies 08 12331 g028
Batenburg et al. [55] developed a semi-steady-state 2D model for UCG in open channels for the developed gasifier only. Unlike other models, they assumed that oxygen instantaneously reacts with the combustible gases present in the channel instead of reacting with the coal surface. Their interest was only to investigate the process within the channel after the injection gas percolates through the inert permeable rubble zone. Reaction rates were calculated based on resistances in the system including boundary layer, pore diffusion, surface phenomena, and chemical kinetics. For laminar and turbulent natural convection along a smooth vertical wall, they used the Sherwood number (Sh) that is correlated with the Grashof number (Gr) and Schmid number (Sc) as follows:
Sh = 0.59 Gr1/4Sc1/4,    Gr < 109   (laminar)
Sh = 0.13 Gr1/3Sc1/3,    Gr > 109   (turbulent)
For forced convection-dominated flow the Sherwood number was expressed in terms of the Reynolds number (Re) as follows:
Sh = 0.027 Re4/5Sc1/3
Heat transfer was modeled by radiation between the walls of the channel. They also included the effect of natural convection due to the temperature difference. Their results showed that the effect of pressure on gas composition is negligible (Figure 29). Their results also indicated that natural and forced convection transfer coefficients are in the same order of magnitude and both of them are important.
Figure 29. Net fluxes and temperature as a function of process pressure for oxygen injection rate = 0.1 mol/m/s and water injection rate = 0.05 mol/m/s [55].
Figure 29. Net fluxes and temperature as a function of process pressure for oxygen injection rate = 0.1 mol/m/s and water injection rate = 0.05 mol/m/s [55].
Energies 08 12331 g029
Pirlot et al. [59] developed a 2D steady-state model by extending the idea of Batenburg et al. [55] with two distinct zones: a low-permeability rubble and ash around the injection point and a high-permeability peripheral zone along the coal wall (Figure 30). After a short transitory starting phase, those two zones were identified by other researchers [64,65] for a thin seam at great depth. During initial combustion, a cavity was identified partially filled with inert materials near the injection hole [64,65]. Their simulation for cavity evolution was based on one main parameter, the permeability ratio between the low- and high-permeability zones. The gasifying agent was assumed to pass through the low-permeability zone surrounding the injection point prior to its arrival in the high-permeability zone, where reactions with the coal wall occur. They combined two separate models: (1) a flow model for the calculation of the flow through the low-permeability porous zone using Darcy’s law and the continuity equation; and (2) a chemical model for the calculation of the chemical processes occurring between gas and the coal wall in the high-permeability zone using empirical correlations for mass and heat transfer for the packed bed by assuming plug flow in the gas phase. The coal consumption rate was calculated only on the channel wall. Their model did not consider the details of moisture and volatile matter released by drying and pyrolysis. They assumed that volatile matter is released in the form of water and hydrogen in proportion to the consumption rate of carbon. They concluded that the permeability ratio is one of the main parameters for determining the success of underground coal gasification because of the observation of the increasing final gasifier area, power, and trial duration with the increase of the permeability ratio. According to them, the cavity growth and shape obtained from their model were in reasonable agreement with the Pricetown field trials.
Figure 30. Conceptual sketch of the two separate zones and species flow direction (redrawn from Pirlot et al. [59]).
Figure 30. Conceptual sketch of the two separate zones and species flow direction (redrawn from Pirlot et al. [59]).
Energies 08 12331 g030
Kuyper [56,60] developed a 2D model to describe UCG process in a cross-section of an open channel (Figure 31) for a typical western European coal layer of thin seams (1–2 m). Field trials of UCG indicated the growth of the cavity upwards and radially outwards around the injection well as gasification proceeds [66]. As a result, at thin seams, the top wall is exposed to rock materials and failure of the overburden is apparently expected. Such a failure would create open channels around the rubble materials deposited on the floor just about the injection hole; as shown in Figure 31. Considering this fact, in their model, the top and bottom walls were assumed as impermeable and adiabatic rock materials. However, the main focus of their work was to obtain an insight for heat and mass transfer due to double-diffusive turbulent natural convection in which both the temperature and concentration gradients play a role in the transport process. The justification of using double-diffusive turbulent natural convection was based on the expectation that Gr/Re2 >> 1, considering the Grashoff number (Gr) of the natural convection flow in the process is approximately 1010 while the Reynolds number (Re) of the forced convection flow is about 104. In this situation, if the forced convection flow is superimposed on the natural convection flow, the heat and mass transport phenomena inside the cavity will be dominated by double-diffusive natural convection and a spiral flow can be expected. The justification for this assumption can be found elsewhere [67]. The Navier-Stokes equation and the k-ε turbulence model were used to describe fluid flow due to large gradients of density and concentration in the channel. Radiation and convection were assumed to be the major heat transfer mechanisms in the channel. Model results showed that oxygen is consumed far from the coal wall by combustible gases. The double-diffusive natural convection in the channel was observed to cause the periodic generation and collapse of CO2 bubbles (Figure 32). Also, mass transfer was reported to be the controlling mechanism (the rate-limiting step) for reduction reactions at the coal wall. They also studied the effect of CO2 injection into the coal seam and reported that CO2 injection has a similar effect as adding steam, although to a lesser extent.
Figure 31. Channel formation in thin seams (Modified from Kuyper et al. [60]).
Figure 31. Channel formation in thin seams (Modified from Kuyper et al. [60]).
Energies 08 12331 g031
Figure 32. Production of CO and CO2 as function of time for κ = 0.12 m−1 and an oxygen injection rate of 47 mmole m−1s−1 [60].
Figure 32. Production of CO and CO2 as function of time for κ = 0.12 m−1 and an oxygen injection rate of 47 mmole m−1s−1 [60].
Energies 08 12331 g032
Because of the identification of UCG application in thick coal seams and the presence of a porous bed of ash overlying the injection point and a void space above the UCG cavity based on field trial excavation [66,68,69] as depicted in Figure 33, Perkins and Sahajwalla [61] considered a thick coal seam and expanded Kuyper’s model by including an ash layer at a lower part of the channel (Figure 34). They developed a 2D axisymmetric model by using Fluent to investigate double-diffusive natural convection along with relevant reactions in a partially filled cavity. Although the bottom wall was considered to be inert and adiabatic, the side and top walls were considered to be carbon due to the thick coal seam. The temperature of the gas and the porous ash bed were assumed to be in thermal equilibrium. The flow though the porous ash bed was modeled using a laminar approximation; however, the turbulence model is only applied to the void region of the cavity. The low Reynolds number turbulence model of Launder and Sharma [70] was employed for the turbulent model by replacing the constants in the standard k-ε model with the function of the turbulent Reynolds number in order to extend the model’s applicability to low-turbulence regions that are close to walls. The following turbulent kinetic energy (kt) and the turbulent dissipation rate (ε) equations were used in order to account for turbulent intensity:
t ( ρ g k t ) + . ( ρ g v ¯ k t ) = . [ ( μ g + μ t σ k t ) k t ] + G k t + G b ρ g ε + Y m + w k t ˙
t ( ρ g ε ) + . ( ρ g v ¯ ε ) = . [ ( μ g + μ t σ ε ) ε ] + C 1 ε ε k t ( G k ) C 2 ε ρ g ε 2 k t + w ε ˙
Figure 33. Schematic of an underground coal gasification cavity from field observation (redrawn from Perkins and Sahajwalla [61]).
Figure 33. Schematic of an underground coal gasification cavity from field observation (redrawn from Perkins and Sahajwalla [61]).
Energies 08 12331 g033
Figure 34. Channel geometry in Perkins and Sahajwalla [61].
Figure 34. Channel geometry in Perkins and Sahajwalla [61].
Energies 08 12331 g034
For other conservation equations for void space, they followed the standard governing equations developed for porous media with the exception of φ (porosity) = 1. For the simplicity of their model, only three chemical reactions (reactions 2, 8, and partial char combustion) were considered. They found that the flow behaviour in the void space is dominated by a single buoyant force due to the temperature gradient resulting from the combustion of oxygen with CO produced from the gasification of CO2 at the walls. From the observable entity, by changing the injection point, they concluded that oxygen should be injected at the bottom of the channel, otherwise valuable gasification products would be oxidized, leading to a low heating value of the production gas. The model was verified with Biezen’s experiments [71] for double-diffusive natural convection in a trapezoidal channel.
Considering the importance of the flame front propagation in the gasification channel for viable operations in UCG, Saulov et al. [72] developed a simplified physical model to describe the primary features of the flame behaviour in a long channel through a coal block. According to their model, a number of factors, such as air flow rate, flame temperature, oxidizer diffusion rate, and radiative heat transport, etc., are the decisive factors to determine if the flame tends to propagate toward the injection point (reverse combustion) or the downstream direction (forward combustion). The flame speed is affected by the air flow rate through the energy balance in the channel. Their model predicts the flame propagation toward the injection point (upstream) for a low air flow rate and toward the downstream direction for a high air flow rate. These predictions from their model are in agreement with qualitative observation by Chernyshv [73]. The flame speed is affected by the air flow rate through the energy balance in the channel. The introduction of oxygen instead of air is found to increase the flame temperature and reactions rates and, as a result, the flame propagates at a faster rate. They also reported a quantitative comparison with the experimental data reported in the book by Skafa [20]. Figure 35 demonstrates two burn velocities, i.e., upstream end velocity, s, and downstream end velocity, sb, as functions of the injection air speed in the upstream undisturbed section of the channel. It is noted that positive values of the velocity mean that the flame front propagates downstream, while negative values imply that the front moves upstream.
Figure 35. Dependence of burn velocity on the air flow rate [72]: (a) experimental data from Skafa [20], and (b) calculated data from Saulov et al. [72].
Figure 35. Dependence of burn velocity on the air flow rate [72]: (a) experimental data from Skafa [20], and (b) calculated data from Saulov et al. [72].
Energies 08 12331 g035
In the case of very hot and diffusion-controlled flame, Saulov et al. [72] considered the limit of high temperatures, a high activation energy, and a strong air flow. Under these conditions the surface of the channel is considered to have two zones, cold and hot. The temperature is found insufficiently high in the cold zone to initiate reactions, while in the hot zone any oxygen on the surface exhausted quickly due to instant reactions. These zones are separated only by a very small distance due to high activation energy. The overall reaction rate is determined by the rate of diffusion of oxygen to the hot zone, while the oxygen concentration on hot walls is essentially nil. Under such conditions the flow in the gasification channel is turbulent and the turbulent flame is fully controlled by diffusion. The injection rate is found to have no control over the flame position. On the other hand, if the combustion of coal begins with devolitalization reactions at low temperatures and these reactions play a noticeable role in initiating the rest of the oxidation process or in the overall energy balance, the flame position is affected by the air speed and becomes controllable. On the other hand, if the devolitalization reactions during the combustion of coal begin at low temperatures and play a noticeable role in initiating the rest of the oxidation process or in the overall energy balance, the flame position is affected by the air speed and becomes controllable. In this case, the flame is not very hot and the cooling effect of the air flow is strong.
The consideration of natural convection is found to be one of the main phenomena in the channel model development. Natural convection plays an important role for the mixing of injected blast gas and the gases coming from the channel wall. The channel model is found to better calculate sweep efficiency. However, most of the channel models neglected drying and pyrolysis which are considered very important in the coal block model. In order to determine cavity shape and size, the channel model is preferred.

2.3. Coal Slab Models

Some models have attempted to describe the UCG coal seam as a coal slab. These models describe the process by movement of various defined regions in a coal slab perpendicular to the flow of the injected blast gas. These regions usually include the gas film, ash layer, char region, dried coal, and virgin coal. The existence of different regions is due to the slow heating rate of the UCG. At a very high heating rate, there is a possibility of the coincidence of a drying front with a combustion front [74]. The general framework of these models is shown in Figure 36. Tsang [74] was the first to use this approach considering the observation of the development of drying, pyrolysis, and gas-char reactions zones around the most permeable path in the coal seam. In addition, the profiles of temperature and saturation as well as the direction of heat and mass transfer exhibited from the pyrolysis experiment (Figure 37) of a small coal block performed by Westmoreland and Forrester [75] inspired him to consider necessary assumptions, such as the existence of a drying front, heat flux, and mass transfer direction, etc., for this approach. In that experiment, a constant heat was applied to the coal surface and career gas was supplied along the length of the cylinder. Some researchers considered the virgin coal region as a wet zone due to the presence of water which is separated from the dry zone by the evaporation/drying front. The evaporation of water is assumed to take place entirely at the retreating drying front. Thus, the model following this approach must be a moving boundary value problem with phase change which is known as the “Stefan” problem. In this approach, as can be speculated from Figure 36, there is an efflux of steam, devolatilized materials, and “self-gasification” products from the wall to the cavity, while there is a counter-current flux of heat towards the cavity wall. “Self-gasification” is considered as the reaction of the gases (steam and devolatilized gases) with hot char while they pass through it before flowing into the cavity. Because of the consideration of different layers, unlike other types of models, for each layer separate mass and energy balances are usually considered. As a result, the governing equations for mass and energy balances are of split boundary types. The mass flux is considered to be diffusion dominant. Table 16 and Table 17 represent a glimpse of a list of coal slab models discussed here with their essential features and reaction rate control mechanisms of the main reactions in the gasification process, respectively.
Figure 36. Scheme of coal block models (modified from Perkins and Sahajwalla [76]).
Figure 36. Scheme of coal block models (modified from Perkins and Sahajwalla [76]).
Energies 08 12331 g036
Figure 37. Schematic diagram of temperature and saturation profile and heat and mass transfer in coal block (redrawn from Tsang [74]).
Figure 37. Schematic diagram of temperature and saturation profile and heat and mass transfer in coal block (redrawn from Tsang [74]).
Energies 08 12331 g037
Table 16. Some essential features of reported coal slab models.
Table 16. Some essential features of reported coal slab models.
ResearchersDimension and Time DependenceHeat TransferMass DiffusionFluid FlowThermo Mechanical FailureWater InfluxHeat Loss
ConductionConvectionRadiation
Tsang [74]1D and T
Massaquoi and Riggs [77,78]1D and S D
Park and Edgar [79]1D and T D
Perkins and Sahajwalla [76,80]1D and PSNS
D = Darcy flow, NS = Navier-Stokes, S = Steady State, T = Transient, PS = Pseudo-steady.
Table 17. Reaction rate control mechanisms of the main reactions in the gasification process.
Table 17. Reaction rate control mechanisms of the main reactions in the gasification process.
ResearchersDryingPyrolysisChar ReactionWater-Gas Shift ReactionGas Phase Reaction
C + O2 → CO2C + CO2 → 2COC + H2O → CO+ H2C + 2H2 → CH4CO + H2O↔ CO2+ H2CO + 0.5O2 → CO2H2 + 0.5O2 → H2OCH4+2O2→ CO2+2H2O
Tsang [74]HP PP K
Massaquoi and Riggs [77,78]HP PP EII
Park and Edgar [79]HPDPP II
Perkins and Sahajwalla [75,80]HP PPPE I
D = Diffusion limited, E = Equilibrium, H = Heat transfer limited, I = Infinite rate, K = Kinetic (power law) and bulk diffusion, P = Power law kinetics.
As stated earlier, Tsang [74] was the pioneer of the approach and developed a 1D unsteady UCG model considering the two regions of a coal block: the wet zone, and dry zone. In the wet zone, heat conduction and liquid transfer are solved. It was assumed that all drying happens in one moving point (Stefan model). In this model, the movement of the drying front was described using the following equation:
 φρ w ( Δ H v ) x ¯ t =   K ¯ d T d x K w T w x ,  T = T vap  at x = x ¯
where ΔHv is the heat of the evaporation of water, x ¯ is the location of the drying front.
This model neglected the effects of steam convection and assumed that the drying process was dominated by conduction. In the wet zone, the liquid diffusion was considered as a function of saturation; however, for the dry zone, multicomponent diffusion was calculated using viscous flow. In the dry region, the pseudo-steady-state assumption was adopted to solve mass balance, heat transfer, and gas flux. Pyrolysis reactions were modeled as simultaneous independent reactions with kinetic parameters obtained from the works of Campbell [49]. Porosity and permeability were set to change with the extent of the reactions with correlations obtained from oil well acidization. This model neglected the combustion and provided the heat by setting a high temperature at the boundary. The experimental results of Forrester [81] were used for model validation. The calculated drying front locations with time were found to be in good agreement with the experimental data (Figure 38). However, the effect of the heating rate on the location of the drying front shows that a higher surface heating rate causes a shallower penetration of the drying front into the coal block (Figure 39). This observation indicates the possibility of the coincidence of the drying front with the combustion front at a very high heating rate. Finally, it was concluded from the model that drying cannot be neglected in the UCG model as it is a factor that dominates the heat transfer rate.
Figure 38. The location of the drying front as a function of time [74].
Figure 38. The location of the drying front as a function of time [74].
Energies 08 12331 g038
Figure 39. Effect of the heating rate on the location of the drying front [74].
Figure 39. Effect of the heating rate on the location of the drying front [74].
Energies 08 12331 g039
Massaquoi and Riggs [77,78] extended Tsang’s 1D model by including the combustion of char and volatiles at the boundary while solving for a steady-state case. Beside the zones shown in Figure 36, they included another zone of bulk gas that is composed of water vapor, CO2, O2, and some inert gases outside the gas film. The developed model was used to describe the simultaneous combustion and drying of a wet Texas lignite coal. The flow of gases in the porous dry coal was modeled by Darcy’s law and heat transfer was considered by both conduction and convection. However, in the wet zone, only conduction was considered for heat transfer. Radiation heat transfer was also considered between the edge of the ash and the bulk gas. Because of the drop in the total pressure across the dry coal and the high permeability of pyrolyzed coal, mass flux was assumed to happen mainly by viscous flow. Multi-component diffusion was assumed to occur at different zones. They considered almost all the reactions listed in Table 1 except for methanation. Following the approach of Tsang [74], pyrolysis reactions were modeled as simultaneous independent reactions with kinetic parameters obtained from the works of Campbell [49]. Because of steady-state assumptions, all the governing equations were written in ordinary differential equations in terms of the moving coordinates. The cavity wall and the drying front movements were considered to have same velocity. In their model, two different cases for combustion were assumed: (i) homogeneous char combustion with the flame located at the ash; and (ii) heterogeneous char combustion with the flame located at the coal surface, and finally a criterion was defined for the transition from heterogeneous to homogenous combustion. For the case of Wyodak subbituminous coal, their result indicated that the criteria principally depend on the conditions of bulk oxygen concentration, bulk gas temperature, coal surface temperature, and coal burning rate. However, different solutions can be obtained at different conditions. For example, as seen in Figure 40, it is found that, for a given bulk gas temperature, there is a limiting bulk oxygen concentration and an oxygen concentration above the limiting bulk; for a fixed oxygen concentration and bulk gas temperature there are two possible solutions (Figure 40): one with the flame at the coal face (at a lower surface temperature) and another with the flame in the ash (at a higher surface temperature). However, heterogeneous combustion was considered unstable due to the decrease of the coal face temperature with the increasing bulk gas temperature for a constant oxygen concentration. In contrast, homogenous combustion was considered stable because of the experimental observation reported in the literature [82,83,84]. They also indicated a maximum value of the bulk gas temperature above which the flame is always away from the coal surface.
Figure 40. Effect of bulk gas temperature of the combustion of wet Wyodak coal [78].
Figure 40. Effect of bulk gas temperature of the combustion of wet Wyodak coal [78].
Energies 08 12331 g040
In addition, they have concluded that the burn rate would be nearly linear with the increase in the concentration of oxygen when the flame front is located in the char face, and the burn rate increases when the flame is located in the gas film (Figure 40). In their model, the effect of the ash layer thickness, coal moisture content, and film transport coefficients were also discussed. A higher film mass transfer coefficient and a lower moisture content and ash layer are observed to shift the transition line towards higher coal surface temperatures corresponding to the higher burning rate for a constant bulk oxygen concentration. For the homogeneous combustion model, the performance of the model was in good agreement with the experimental data [83,84].
Park and Edgar [79] developed an unsteady 1D model with a moving burning front based on the work of Massaquoi and Riggs [77,78] to describe lateral cavity growth in UCG. However, unlike the assumption of Massaquoi and Riggs [77,78], they did not consider having the same velocity for the cavity wall and drying front during the early stage. Beside char gasification, they included coal shrinkage due to drying and pyrolysis, as well as the “self-gasification” of drying and pyrolyzed products (steam and CO2) in the region between the cavity wall and the drying front in order to consider the additional movement of the cavity wall. The cavity wall movement was determined by the rate of the removal of coal by chemical reactions and the physical movement of the cavity wall due to the shrinkage of the coal. In their simulation, an increase of the cavity growth was noticed during the transient period due to the shrinkage of the coal; however, the cavity growth with and without considering the shrinkage eventually merged into one rate when the steady state was reached (Figure 41). However, they suggested that the movement of the cavity wall due to shrinkage is only important in lab-scale processes and can be neglected in larger-scale processes because of the long timescales in which the transient period is negligible compared to the steady-state period. Their results indicate that cavity growth is controlled by the rate of oxygen transfer to the cavity wall, when the flame is located at the char surface. If the flame is separated from the cavity wall, CO2 and steam gasification determine the cavity growth rate. The moisture content in coal was found to be a complex parameter for cavity growth. Although the surface temperature decreased with increasing moisture content, the higher moisture content indicated a higher cavity growth at a high coal surface temperature (above 1000 °C). This observation implies that the availability of water in the coal, rather than the reaction temperature, controls the steam gasification reaction rate. They validated their results with the experimental data of lateral cavity growth from a Texas lignite core combustion performed by Poon et al. [85]. From the models discussed above [77,78,79], it can be speculated that under oxidizing conditions, the cavity growth rates are controlled by mass transfer.
Figure 41. Effect of the shrinkage on the cavity growth rate [79].
Figure 41. Effect of the shrinkage on the cavity growth rate [79].
Energies 08 12331 g041
Perkins and Sahajwalla [76,80] also developed a 1D transient coal block model with an extension of Tsang’s study [74] by including multicomponent diffusion and the random pore model to account for changes of heterogeneous reaction rates with conversion. For multicomponent diffusion, Stefan-Maxwell equations and the bi-dispersed dusty gas model were used for the gas film and dried coal matrix, respectively. Because of experimental evidence (Figure 42) of different sizes of pore distribution in the porous structure of coal, a bimodal porous structure (macropores and micropores) was considered where the macropores act as the major conduits for species flow and the micropores act as the dominant region for gasification reactions.
Figure 42. Pore structure of raw American coals on dry and ash-free basis [70].
Figure 42. Pore structure of raw American coals on dry and ash-free basis [70].
Energies 08 12331 g042
They have proposed that cavity growth occurs at reducing conditions, therefore only heterogeneous gasification reactions are solved, and required heat was provided by defining a constant temperature at the char surface. The char surface was allowed to participate in radiation heat exchange with its surroundings. For pyrolysis, they followed the work of Tsang [74]. The movements of the drying front and char front were assumed to be equal under pseudo-steady-state conditions. They assumed that solid and gas phases are in thermal equilibrium and bulk gas has a fixed composition that is representative of the product gas. In their first paper [80], they validated their model with the experimental data obtained by Forrester [81]. The drying front location as a function of time was in good agreement after using a scaling factor of 1.2 for thermal conductivity (Figure 43). The temperature profile through the block at different surface temperatures also agreed well with the experimental data. However, when the cavity growth rates as a function of gas temperature for two different bulk gas compositions were compared, no significant difference was observed (Figure 44).
Figure 43. Drying front location as a function of time for cylindrical coal block heating experiments [80].
Figure 43. Drying front location as a function of time for cylindrical coal block heating experiments [80].
Energies 08 12331 g043
Figure 44. Cavity growth rate, as a function of gas temperature, for two different bulk gas compositions [80].
Figure 44. Cavity growth rate, as a function of gas temperature, for two different bulk gas compositions [80].
Energies 08 12331 g044
In a later paper [76], they extended their work by including the ash layer and water influx and investigated the effect of various parameters and operating conditions (e.g., temperature, pressure, water influx, gas composition) on the process performance and cavity growth rate. In order to compare the relative effect on cavity growth due to operating conditions and coal properties, they explored the sensitivity of the cavity growth rate (vcv) to changes in parameters (pr) using the following finite difference formula:
s ( p r ) = v c v ( p r o + Δ p r )   v c v ( p r o Δ p r ) 2 v c v ( p r o ) p r o Δ p r
where s(pr) is the sensitivity of the parameter, pr° is the base value of the parameter, and Δpr is a perturbation to the value of the parameter.
A ranking was established according to the sensitivity of the parameters to the cavity growth rate. A positive value was considered as an increase in vcv, while a negative value was considered in reverse order. Table 18 and Table 19 show the sensitivity of the cavity growth rate to perturbation in operating conditions and coal properties around the base values, p°. Since the summation of the composition of proximate analysis is equal to one, changing any one of these parameters was followed by adjusting the other parameters in order to maintain the constraint of the summation of one.
Table 18. Sensitivity of the cavity growth rate to changes in operating conditions [76].
Table 18. Sensitivity of the cavity growth rate to changes in operating conditions [76].
ParameterΔp/p° = 10−3Δp/p° = 10−2Rank
Bulk gas temperature, K1.5 × 1035.2765.2761
Water mass flux, kg/m2·s1.0 × 10−3−0.114−0.1152
Mole fraction of steam in bulk gas2.0 × 10−10.0660.0713
Ambient pressure, Pa1.5 × 1060.0790.0674
Mole fraction of CO2 in bulk gas1.2 × 10−1−0.009−0.0095
Mole fraction of H2 in bulk gas1.6 × 10−10.0080.0086
Heating rate, K/s1.0 × 10−20.0060.0077
Mole fraction of CO in bulk gas8.0 × 10−20.0050.0058
Table 19. Sensitivity of the cavity growth rate to changes in coal properties [76].
Table 19. Sensitivity of the cavity growth rate to changes in coal properties [76].
ParameterΔp/p° = 10−3Δp/p° = 10−2Rank
Critical conversion of carbon at which the coal spalls9.50 × 10−1−3.719−3.7411
Mass fraction of fixed carbon from proximate analysis3.35 × 10−1−1.359−1.3702
Solid material density (kg/m3)1.47 × 103−0.913−0.9253
Mass fraction of ash from proximate analysis2.40 × 10−10.7140.7284
Ash layer thickness (m)5.00 × 10−3−0.689−0.7035
Mass fraction of volatile matter from proximate analysis3.00 × 10−10.5890.5766
Mass fraction of moisture from proximate analysis1.25 × 10−10.3440.3317
Emissivity of ash (for radiation)8.50 × 10−10.1140.1278
Initial specific surface area (m2/m3)1.23 × 1070.1050.1059
Results of their parametric study indicate that temperature of the coal surface, water influx, and coal composition have the highest impact on the cavity growth rate. Decreasing fixed carbon was found to increase the cavity growth rate. This implies that the cavity growth rate is more for low-rank coal due to the smaller fraction of fixed carbon in the proximate analysis. Perkins and Sahajwalla [76] graphically demonstrated the variation of cavity growth due to the variation of different parameters through their modeling. In addition, they incorporated the coal spalling model by changing the value of the critical conversion, Xcritical. Due to lack of information on coal spalling, they considered the critical conversion of fixed carbon as an indicator of coal spalling by assuming a level of carbon conversion at which some fixed carbon and volatile matter are unconverted and detach from the coal block due to mechanical stress and breakage. Critical conversion was found sensitive to the bulk gas temperature, molar flux of product gas, and water influx.
Finally, they compared the estimated cavity growth rate from their simulation with six different field trials (Hanna I, Pricetown, Large Block Test No. 5 (LBK-5), Partial Seam CRIP (PSC), Rocky Mountain I (RM I), and El Tremedal) using the available operating conditions and coal properties of those field trials and obtained a good agreement with the data from the field trials (Figure 45) despite some uncertainty of the field data. Although they conducted a detailed study of the coal block model, they did not consider char combustion and other gas phase reactions (except the water-gas shift reaction) which are strong functions of the temperature. Rather, the assumption of a constant temperature and fixed gas composition at the boundary limits the applicability of the model to field trials.
Figure 45. Comparison of estimated cavity growth rates from six UCG field trials with model simulations [80].
Figure 45. Comparison of estimated cavity growth rates from six UCG field trials with model simulations [80].
Energies 08 12331 g045
The special feature of the coal slab model is the tracking of the drying and combustion front movement. This model can successfully demonstrate the drying and devolatilization behaviour of large coal particles; however, this model is yet to be validated using UCG field trial data. Due to the assumption of semi-infinite coal slab, it can be speculate that this model is only good for a thick coal layer. All the models developed thus far by considering the coal slab are one-dimensional (1D). As a result, information on cavity shape cannot be obtained from these models.

2.4. Other Approaches for UCG Modeling

2.4.1. Reactor Models

Over the years, there have been a few studies in which the UCG process is described as a series of ideal reactors. Such models can be used to predict product gas composition without considering complex phenomena occurring underground.
Chang et al. [86] used this approach to predict the resource recovery during the UCG process. In this model, the domain was divided into cylindrical elements where each element consists of two distinct zones: rubble and void. It was assumed that the rubble zone acts as an ideal plug flow reactor (PFR), while the void space was assumed to be perfectly mixed as a continuous stirred tank reactor (CSTR). Each element was allowed to change based on coal consumption. The fraction of flow in each zone was changed to match the results of the field trials. The spalling and water influx were calculated based on empirical correlations. Their modeling results highlight the importance of turbulence in field trials. Therefore, they concluded that small laboratory-scale experiments could not be scaled to field trials.
In order to obtain in situ measurements, the stimulus response technique (tracer injection) was used in some field trials to study the flow distribution inside the underground cavity. Some researchers have used the response of tracer material from the experiment to characterize different zones in the underground cavity and track the cavity growth. Based on the residence time of tracer elements, the cavity is modeled as an arrangement of ideal reactors.
Debelle et al. [87] used various chemical engineering reaction models to fit the quasi-exponential decay of the tracer response from a UCG experiment in Thulin, Belgium. In the experiment, 133Xe was selected as a tracer because of its detectability from a complex gaseous mixture by using a γ-ray detector with high precision. The quasi-exponential decay of the tracer was considered as an indication of mass exchange between a flowing fluid and a dead zone. That is why the models were developed to simulate the flow conditions based on the material exchange between flowing fluid and the dead zone. The film theory was applied to simulate the mass transfer between a stagnant volume fraction and the following fluid volume fraction. However, for the interaction between the flowing fluid and porous zone, molecular diffusion was considered. The dead zone was represented either by a stagnant zone or by a homogenous porous zone where materials exchange with the flowing fluid. Various stages of the process such as reverse combustion linking and forward gasification were represented using regression models. For the gasification period, it was concluded that the cavity offers a high level of back-mixing, resembling a well-mixed reactor with a small dead zone.
Thorsness [88] performed tracer experiments in a UCG field test at the Hoe Creek trial III. Their goal was to investigate the extent of the water influx and steam injection requirements. The tracer response curve was regressed with a series of CSTRs in two parallel paths.
Pirard et al. [89] considered the El Tremedal trial in Spain where helium was used as a tracer. In their study, two classical engineering models were used: (i) the combined model composed of a succession of stirred tanks in a series exchanging matter with an adjacent porous zone (the STEP model) and (ii) the combined model composed of a succession of stirred tanks in a series exchanging matter with a stagnant zone (STES model). Both models are three-parameter models. The parameters are: the volume fraction occupied by the flowing fluid, the number of tanks in the series, and the dimensionless parameter characteristic of the pore diffusion. By these methods, they were also able to detect cracks or breaks in the coal seam. Their results indicated the STEP model as the best fitted model, and finally they concluded that the cavity behaves almost like a small number of CSTRs in a series with a high level of back-mixing.
Recently, Daggupati et al. [90] used computational fluid dynamics (CFD) to describe the cavity as series of ideal reactors (Figure 46). The residence time distribution (RTD) of gases was calculated based on CFD results. This approach was validated by laboratory-scale tracer experiments. Based on RTD and velocity distribution, the cavity was simulated as a network of ideal reactors. The proposed arrangement of reactors consists of a series of five CSTRs and a side stream of five small CSTRs, followed by a plug flow reactor (PFR). They concluded that as the cavity size increases, the cavity behaviour changes from PFR to CSTR. They have implemented their reactor model in ASPEN PLUS (chemical process optimization software) and compared the results with CFD predictions.
Figure 46. Proposed reactor arrangement for the cavity in the works of Daggupati et al. [90].
Figure 46. Proposed reactor arrangement for the cavity in the works of Daggupati et al. [90].
Energies 08 12331 g046

2.4.2. Probabilistic Simulation

In a unique approach to model the UCG process, Batenburg and Bruining [91] adapted a probabilistic simulation to describe the forward gasification phase in a 2D UCG geometry (Figure 47). The model used stream function on the grid blocks. The derivative of the stream function along the ash-char interface was interpreted as the possibility that the interface was moving in a certain block. This model is very limited, as the temperature of the gases and composition are constant in each region.
Figure 47. Probabilistic simulation of the UCG process by possibility of movement of the following interfaces: 1.void-rock; 2. coal-rock; 3. void-coal; 4.void-ash; 5. ash-void [92].
Figure 47. Probabilistic simulation of the UCG process by possibility of movement of the following interfaces: 1.void-rock; 2. coal-rock; 3. void-coal; 4.void-ash; 5. ash-void [92].
Energies 08 12331 g047
From the same research group, Biezen et al. [93] extended the concept of the probabilistic simulation by including the movement of several interfaces including ash-void, ash-coal, void-coal, and void-rock. The model consists of two modules: one module solves flow in the entire domain and the other selects a block of coal/rock for gasification or roof spalling. However, in order to consider roof spalling, they incorporated a very simple spalling model which assumes a constant spall rate. A single streamline was selected randomly and based on the interface located on the streamline, and relevant physical phenomena proceeds. For example, in the void-coal interface surface gasification reactions occur and the corresponding amount of ash accumulates at the bottom of the cavity. They later extended their model to simulate cavity evolution in a 3D domain [92]. This model was intended to simulate the development of a gasifier from early times (gasification in a permeable bed) to a developed stage (gasification in a channel). They did not solve the mass and energy balance together. As a result, the model is limited for making predictions of performance. However, the model is more capable of interpreting previous field trials. The application of this model is limited due to its complexity and unphysical assumptions such as constant temperature at each domain.

2.4.3. CAVSIM (a computer code for 3D, axisymmetric UCG model) Model

All of these reviewed modeling efforts so far have the deficiency of focusing on some aspect of the UCG process while neglecting others. During UCG trials in the US, LLNL funded many modeling studies for the UCG process. Several 1D models were developed and validated with certain series of data from field trials. These segregated models were later combined together to form the CAVSIM process model that represents 15 years of continuous UCG research and development in the US [23]. The model is applicable for predicting lateral cavity growth of thick and shrinking coal seams in which oxygen is injected at the bottom of the coal seam. The model was based on some major assumptions such as: the cavity was axisymmetric around injection; thermal radiation was the main heat transfer mechanism in a well-mixed void space; cavity growth was dominated by the thermo-mechanical failure of the wall, and a packed bed of char and rubble forms over a thin layer of ash. Various sub-models were incorporated into this model to quantify major phenomena in UCG such as water influx, flow dispersion through a rubble bed, sidewall growth of the cavity due to reactions, and the spalling of the cavity roof into the domain. This model reproduced results of some of the field trials in the US by tuning several parameters that were introduced in the model. Therefore, the model could not be used for predicting UCG behaviour in other seams, as these parameters are not known prior to the process.

2.5. CFD Software Tools

From the above reviews, it is evident that the understanding of gasification for a particular coal seam in a specified location with its boundary conditions and various time-varying phenomena are imperative in order to identify the appropriate UCG model. Although several models with varying levels of complexity have been published over the past years, the applicability of these models is limited to specific and isolated cases and the models could hardly be used to predict the performance of UCG trials. Furthermore, most of these models are in 1D or 2D, while the field trials reveal a 3D non-regular cavity shape. Thus, it is of the utmost importance to develop a 3D model in order to predict the performance of UCG trials. The earlier models compromised between model complexity and simplicity in order to get faster computation runtimes. However, currently, with the advent and advancement of different modeling software tools and powerful computers, it is possible to develop a 3D model with a varying limit of complexities and a much faster runtime. There are several commercial CFD software tools which were used successfully for the UCG models. Some of the popular tools are:
-
Ansys Fluent.
-
CMG Stars.
-
Star CCM+.
-
Comsol Multiphysics
-
ABAQUS
The basic difference of the above tools is the discretization approaches. Ansys Fluent and Star CCM+ are finite volume-based solvers whereas CMG Stars is a finite difference-based solver. Comsol Multiphysics and ABAQUS follow the finite element approach. Recently, a number of investigators used CFD software tools in order to develop the UCG models. Table 20 shows a list of the researchers who developed a UCG model using different software tools.
Table 20. List of the researchers who used CFD software tools in order to develop their UCG model, recently.
Table 20. List of the researchers who used CFD software tools in order to develop their UCG model, recently.
CFD Software ToolsInvestigators
3D Model2D Model
Ansys FluentShirazi et al. [50]Perkins and Sahajwalla [61]
Luo et al. [58]
CMG StarsKariznovi et al. [9]
Nourozieh et al. [12]
Seifi et al. [94]
Star CCMNitao et al. [95]
Comsol Multiphysics Shirazi [27]
Srikantiah et al. [96]
ABAQUSYang et al. [97]
Sheng [98]
Sarhosis et al. [99]
The commercial CFD software packages are intended for general use. As a result, not all the functions required for the development of UCG models are available in those packages. However, most of the packages have the facility to incorporate user-defined functions. Investigators usually take this benefit and develop their own user-defined functions for their model in order to overcome this shortcoming of the software packages.

2.5.1. 3D Models Using CFD Software Package

There are few 3D models developed so far in the context of UCG. As mentioned earlier, Biezen [92] developed a 3D model for cavity growth in 1995 that appears to be the first 3D UCG model. However, his model was limited due to the absence of any heat transfer calculation and the assumption of a constant reactor temperature. Recently, due to the rapid development of commercial simulation tools for handling complex transport phenomena and chemical kinetics, the development of a 3D UCG model using commercial software tools is especially noticeable. Nourozieh et al. [12] developed a 3D model for CRIP configuration using CMG Stars, a reservoir simulation software, and discretized the domain into coarse cubic meshes of 25–50 cm, typical to reservoir solvers. Their computational domain was similar to the proposed field trials in Alberta (15 m × 5 m × 9 m). They have introduced an arbitrary clay layer in the seam to investigate the effect of inert layers on the process; however, thermo-mechanical failure of this layer is not included in the model. They included chemical reactions (all the reactions listed in Table 1 except reaction 9) and the conduction/convection of heat and species in their model. They did not make any attempt to treat heterogeneous reactions separately. Except for the methanation reaction, all other reactions were modeled using the power law and the temperature-dependent rates were described by the Arrhenius correlation. Their results show a large degree of grid dependence. Resizing grids affected the temperature distribution and altered the gas composition at the producers (Figure 48 and Table 21).
Figure 48. Effect of the grid size on the temperature distribution: large grids (left) and small grids (right) (reproduced from Nourozieh et al. [12]).
Figure 48. Effect of the grid size on the temperature distribution: large grids (left) and small grids (right) (reproduced from Nourozieh et al. [12]).
Energies 08 12331 g048
Table 21. Effect of the grid size on product gas composition [12].
Table 21. Effect of the grid size on product gas composition [12].
Gas Composition (Mole Fraction)CO2H2COCH4
Grid size0.5 × 0.5 × 0.5 m0.34900.28870.22060.1416
0.25 × 0.25 × 0.25 m0.28470.25570.28640.1732
Because of the 3D model, they were able to exhibit the cavity growth and shape (combustion front) that was nearly the shape of tear drops (Figure 49). Their study also showed an increase of gas temperature with an increasing oxygen flow rate (keeping the H2O/O2 ratio constant) and decreasing H2O/O2 ratio (keeping the O2 flow rate constant). With the increasing pressure and H2O/O2 ratio, a significant increase of methane and decrease of hydrogen is noticed in their model prediction. This phenomenon was explained as the favourable condition for methanation at high pressure. However, the methane concentration was found higher even in low pressure conditions as compared to other studies. As they did not consider reaction 9, which converted H2 to H2O and reduced the concentration of H2, the presence of a higher H2 concentration possibly favoured the methanation reaction. They did not present any validation for their model. They have extended their modeling to simulate CRIP with successive ignition points [94]. Again, no validation or comparison with experimental data was reported.
Figure 49. Cavity shape corresponding to the geological structure after 10 days of simulation: x-z cross (right) and x-y cross (left) (reproduced from Nourozieh et al. [12]).
Figure 49. Cavity shape corresponding to the geological structure after 10 days of simulation: x-z cross (right) and x-y cross (left) (reproduced from Nourozieh et al. [12]).
Energies 08 12331 g049
Recently, Shirazi et al. [50] developed a small-scale 3D packed bed model (3 cm × 1.5 cm × 2 cm) for a shrinking coal seam using Ansys Fluent. The following set of governing equations is solved in his model:
Continuity equation:
t   ( φ ρ g ) + ( μ ρ g v ¯ ) =   S
Momentum equation:
t   ( φ ρ g v ¯ ) + . ( φ ρ g v ¯   v ¯ ) = μ 2 φ v ¯ φ p +   φ ρ g g ¯   ( φ 2 μ α v ¯ + φ 3 C 2 2   ρ | v | v ¯ )
The flow in porous medium is implemented in Ansys Fluent by adding a source term into the original Navier-Stokes equation.
Energy equation:
t   ( ρ g C p g φ +   ρ s C p s ( 1 φ ) ) +   ρ g C p g ( v ¯ . T ) = . [ k e f f T ] + Σ r a t e i Δ H i
In the model, it was assumed that gas and solids that form the porous media are in thermal equilibrium in each cell; thus, one energy equation was used to describe heat transport in each cell. Only heat conduction was considered.
Transport of species:
t   ( ρ g Y i φ ) +   ρ g ( v ¯ . Y i ) = . [ ρ g D i Y i ] + Σ r a t e i
Transfer of species was described by the convection-diffusion equation in porous media.
The cavity development is tracked by the decrease of porosity due to chemical reactions and thermal effects. Porosity increases as the solids are being consumed by reactions or species transfer from the solid phase to gas phase due to thermal effects. The change of porosity with chemical reactions is calculated with the following equation:
t   ( φ ) =   S φ
where Sφ is the source term for porosity based on reactions in the coal.
Coal block is considered to be a porous medium with a defined initial porosity and permeability. However, the porous media model in Fluent does not allow for solids to participate in reactions. That is why no separated mass conservation equation was solved for the solid phase. Therefore, in order to simulate and track the solid and gas/solid reactions, a virtual solid phase was defined and a number of user-defined scalars (UDS) were written to be implemented in Fluent. Three scalar equations calculate the consumption of moisture, volatile matter, and fixed carbon following a general form of UDS equations in Fluent:
t ( ρ g φ ) + x i ( ρ g v i x i ) = S
The ash content of the coal is assumed to segregate after consumption of the coal and no ash layer remains over the coal. This assumption has been verified by several experiments and field trials [34,68]. In the model, the changes in porosity and permeability with time due to the drying, pyrolysis, and reactions with coal were included. The gas properties are taken from the Fluent database. Heat capacities of gases were taken following the precise polynomial model.
Their model revealed the effect of different parameters on cavity growth as well as cavity shape and size. Lower oxygen flow rates at the inlet led to a more bulbous cavity shape, whereas at higher flow rates, the cavity became elongated in the direction of the production well. In addition, increasing the inlet oxygen flow beyond a certain value decreased the concentration of CO and H2 in syngas due to the excessive combustion of these gases. Initial permeability of the seam also had a strong influence on the growth rate and shape of the cavity. The higher coal permeability led to faster growth of the cavity and resulted in a much wider cavity (Figure 50).
Figure 50. Cavity shape for cases with different initial permeabilities: (a) 1 mD, (b) 100 mD, and (c) 10 D (reproduced from Shirazi et al. [50]).
Figure 50. Cavity shape for cases with different initial permeabilities: (a) 1 mD, (b) 100 mD, and (c) 10 D (reproduced from Shirazi et al. [50]).
Energies 08 12331 g050
Unlike some model predictions (i.e., Daggupati et al. [90]), a very negligent effect of the temperature of the inlet oxygen was observed on cavity growth. Their model was validated with the lab-scale experiments of Daggupati et al. [34] in terms of cavity growth (Table 22). The agreement between the model calculations and the lab experiments indicates that the relevant physics and reactions of their model for small coal seam geometry are useful in the design of the gasification process.
Table 22. Comparison of the model-predicted cavity dimension with the experiments [47].
Table 22. Comparison of the model-predicted cavity dimension with the experiments [47].
DimensionModel [42]Experiments [25]
Forward length11
Backward length0.4290.433
Height0.3990.517
Width0.4040.425

2.6. Sub Model: Drying and Pyrolysis

The drying and devolatilization are the integral parts of the UCG process. During UCG, drying, pyrolysis, and combustion of the coal occur simultaneously [77,78]. The drying is related to moisture release whereas devolatilization is related to the production of char by the evolution of volatile matters and tars. The development of drying and devolatilization sub-model is also as important as the development of the cavity growth model as all these processes are involved in UCG cavity expansion.

2.6.1. Drying Sub-Model

Drying is usually considered as the release of the unbound moisture in the coal at the evaporation temperature of water. Moisture content is an integral part of the coal seam. However, the percentage of moisture present in coal depends on coal ranks and sources. Although drying was considered with importance in most of the packed bed and coal slab models, this phenomenon has been neglected in most of the channel models. However, instead of a drying rate, the drying front movement calculation is found to be the main focus for coal slab models.
In order to interpret drying phenomena, the following models are well recognized [74]:
-
The diffusion theory,
-
The capillary flow theory,
-
The evaporation-condensation mechanism,
-
Simultaneous energy and moisture transfer in porous media
The evaporation-condensation mechanism is generally found to be used for wet coal where mobile water is available or in a situation where water influx is considered [47,48,100]. In construct, drying is usually considered for the evaporation of moisture content (fixed water) in the coal. This moisture can be present in the coal’s structure as unbound, weakly bound, and chemically bound. The chemically bound moisture is released at higher temperatures around 200–400 °C [101], while the unbound moisture is released at the evaporation temperature of water [102] which is considered as the standard drying temperature for most of the UCG models. Some UCG models assume that drying occurs at a front layer in which the temperature is kept constant, equal to the water evaporation temperature, and drying is assumed to be governed by an Arrhenius-type reaction expression:
d M d t = r a t e d r y i n g = k ( M * M )
where M is the moisture content of coal at any moment and M* is the moisture content of coal from proximate analysis.
However, the release of weakly bound and chemically bound moistures at higher temperatures is not correctly governed in these models. On other hand, some models recognize that the moisture in the coal is distributed throughout the coal and is not located solely at the gas/solid interface. At high pressure, the drying temperature depends on the water vapor pressure. Considering this fact, Thorsness et al. [48] developed an expression for the drying rate that is less than the evaporation rate for mobile water. The evaporation/condensation rate was considered by a pure mass transfer controlled expression as follows [48]:
r a t e e v a p o r a t i o n / c o n d e n s a t i o n = k y ( p w P y w ) ,   for   p w P
where ky is the mass transfer coefficient, yw is the mole fraction of water in gas, pw is the vapor pressure of water, and P is the total pressure. The temperature and vapor pressure are given by:
p w = 10 5 e x p ( 12.61   4690 T d )
T d = β   T s + ( 1 β ) T g
where β is the weighing coefficient. However, most of the models considered β = 1 which reduced the drying temperature to the surface temperature (Td = Ts). From Equation (44), Perkins and Shajwala [100] deduced an expression for the drying temperature for calculating the drying temperature at high pressure:
T d =   4690 12.61 ln ( p w 101325 )
For the drying rate, Thorsness et al. [48] expressed the basic rate for spherical symmetric coal particles as follows:
d M d t r a t e d r y i n g = k y γ 2 3 η ( p w P ρ w ρ w o y w )
where ρ°w is the uniform initial moisture density and parameters γ and η deal with the internal resistance. Parameter η is determined from the following expressions:
η   =   k y ρ w 2 P ρ w o D e f f
Parameter γ is the first Eigenvalue given by:
t a n γ n = γ n ( 1 η )
However, for simplicity, Winslow [47] and Thorsness and Rozsa [42] calculated the drying rate from the evaporation rate described by Equation (43) multiplied by a factor less than one (<1). Due to the inherent dependency of the mass transfer coefficient on particle diameters in the Thorsness [48] model, the drying rate expression also closely depends on the initial size and shape of the coal particles.

2.6.2. Pyrolysis Sub-Model

During devolatilization, the releasing of the volatiles occurs in the following order: chemically bound H2O, CO2, CO, CH4, tars, and H2 [15]. Some other hydrocarbons also evolve at a lower concentration. Due to the high residence time of gas and high pressure in the underground coal gasification process, tar would likely undergo secondary reactions such as cracking, and tar yield would be lower compared to the surface gasifiers [80]. Thorsness et al. [48] used an empirical formula C9Hc for tar and a trace amount of hydrocarbon not accounted for in the volatiles mentioned above. Pyrolysis of coal can be represented by the following reactions:
D r i e d   C o a l C h a r + T a r + A s h + Volatiles
The volatiles can be further represented by:
V o l a t i e l s a 1 C O 2 + a 2 C O + a 3 C H 4 + a 4 H 2 O + a 5 H 2 + a 5   O t h e r   h y d r o c a r b o n s
The matrix of the coefficient can be obtained from the element balances of the atoms involved in the above equation and the data from the proximate and ultimate analysis of the coal. Some models obtained the stoichiometric coefficients based on the work of Cambell [49] or Suuberg et al. [103] for subbituminous coal or lignite coal, respectively.
In order to account for tar, Shirazi et al. [42] followed the kinetic model proposed by Boroson and Howard [90], who considered tar cracking with the following overall reaction:
t a r b 1 C O 2 + b 2 C O + b 3 C H 4 + b 5 H 2
The stoichiometric coefficients and kinetic parameters can be found from the work of Boroson and Howard [104].
The rate of gas evolution is different at different temperatures. Figure 51 shows an experimental observation of gas release from coal as a function of temperature. A similar observation was reported by Campbell [49]. It is observed that different species evolve at different peak pyrolysis temperatures.
Figure 51. Effect of temperature of gas evolution during the pyrolysis of Wyodak coal at a heating rate of 3.33 °C/min under atmospheric pressure on a dry basis (reproduced from Juntgen and van Heek [105]).
Figure 51. Effect of temperature of gas evolution during the pyrolysis of Wyodak coal at a heating rate of 3.33 °C/min under atmospheric pressure on a dry basis (reproduced from Juntgen and van Heek [105]).
Energies 08 12331 g051
There are multiple steps of decomposition during the pyrolysis process. However, in order to avoid complications, Van Krevelen et al. [106] proposed the whole pyrolysis process as a single-step decomposition as follows:
d V M d t =   r a t e p y r o l y s i s = k   ( V M *   V M )
where VM is volatile matter content at any moment and VM* is the effective volatile matter content of coal. VM* is usually assumed to be the same as volatile matter from the proximate analysis. The rate constant k is typically related with the temperature by an Arrhenius expression:
k = k 0   exp (   E R T )
In order to avoid complicacy, some models used the above equation for pyrolysis [9,12,50]. Kariznovi et al. [9] summarized the experimental kinetic parameters for different coal. However, the lumping of a number of independent first-order reactions into a single first-order reaction tends to underestimate the frequency factor, ko, and the activation energy, E. In order to avoid this complexity, Juntgen and Van Heek [105] proposed simultaneous independent reactions with different reaction constants to describe the evolution of different species during pyrolysis as follows:
d V M i d t =   r a t e p y r o l y s i s = k j ( V M j * V M j ) ,   j = 1, 2 .... number of volatiles
where,  k j = k j   exp ( E j R T )
Campbell [49] and Suuberg et al. [103] carried out excellent work to estimate the frequency factor and activation energy for each species for the pyrolysis of Montana lignite and Wyodak subbituminous coal, respectively. A number of researchers [74,76,77,80] included simultaneous independent reactions for pyrolysis with the quantities of pyrolysis gases evolved and the corresponding kinetic parameters obtained from the works of Suuberg et al. [103] or Campbell [49] depending upon the types of coal they used in their models.
Similar to the UCG model, the location of the drying and devolatilization zone depends on the approach of the modeling. Figure 52 shows the reaction zones in the UCG process for different stages of gasification. Usually, a UCG model based on the packed bed approach followed Figure 52a for drying and devolatilization. On the other hand, Figure 52b was followed by the model based on the channel approach.
Figure 52. (a) Reaction zones along the direction of gas flow; (b) Reactions zones surrounding a channel (reproduced from Tsang [74]).
Figure 52. (a) Reaction zones along the direction of gas flow; (b) Reactions zones surrounding a channel (reproduced from Tsang [74]).
Energies 08 12331 g052

2.7. Consideration of UCG Modeling Parameters

As mentioned earlier, UCG is site-specific. Controlling UCG is less flexible as compared to surface gasifiers because some important variables such as the thickness, quality of the seam, moisture contents, and other physical and chemical parameters are unique properties of the natural location and cannot be changed. Most of the models indicated the sensitivity of the process to the properties of the seam and the change of these properties with temperature. In order to obtain the best possible results from governing equations, appropriate solid and gas properties are required. Most of the models discussed above used either constant properties or transient properties obtained from the approximate model or experimental results of others. However, for a realistic model of any specific UCG site, the physical properties of coal and their changes with temperature and pressure, such as thermal conductivity, specific heat, density, porosity, and permeability, etc., need to be considered. Similarly, properties of gas species and their changes with temperature and pressure, such as density, viscosity, thermal conductivity, diffusivity, and heat capacity, etc., need to be considered also. Some experiments might be required in order to formulate these dependencies for the specific seams. In order to enhance the features of the 3D model and to apply the model in UCG field trials, it is imperative to incorporate the coal properties, the change of properties with temperature, and the coal seam geology of the specific site into the model.
In order to avoid model complications and to reduce computation time, some models assumed that gas and solids in the porous media are in thermal equilibrium and that is why they used one combined heat transfer equation. However, there is a large difference between the solid and gas characteristic timescale. It is assumed that gas reaches the steady state before any change happens in the solid [46]. Thus, two separate heat transfer equations for gas and solid, respectively, seem to be more appropriate to incorporate into the model. In some models, the mass conservation of the solid is avoided since drying, pyrolysis, and chemical reactions are considered in mass conservation of gas. However, in order to account for the mass loss of solid as well as ash production, the solid species balance seems to be reasonable for incorporating into the model. The choice of the appropriate use of reaction rates is always challenging due to the complex situation of UCG operations. As these values are obtained from experiments, a number of experiments were performed in order to obtain reaction rates. However, obtaining a reaction rate for a specific reaction needs to consider the operating pressure and temperature similar to UCG operational conditions.
There is hardly a UCG site without any water influx from the surrounding strata. Thus, in a UCG site, water influx is unavoidable and the gasification is performed above the hydrostatics pressure in order to control the water influx; this allowed the water to be used for steam gasification in order to utilize the heat of combustion and increase the concentration of CO and H2 in the product gas. It is imperative to include a water influx formulation that is practicable for a specific site’s coal hydrology. Although some investigators claimed a steam/oxygen ratio of 2.5 as optimum for UCG [34,107], for different sites and coal seams, it is also necessary to explore the appropriate steam/oxygen ratio.

3. Challenges of UCG Commercialization

In spite of the significant advantages of UCG technology and repeated trial successes, UCG is still not close to full commercialization. Historical evidence shows that the commercialization of a technology requires a reasonable balance among the technological, environmental, social, and economic factors. In terms of technology, a lot of problems associated with UCG were solved already in the last few decades. The improvement was focused on obtaining (i) good and nearly constant gas quality and heating value; (ii) high thermal efficiency and resource recovery; (iii) high process efficiency and equipment reliability; (iv) high predictability; (v) good control on the combustion front; (vi) handling bituminous coal; and (vii) handling the depth of the coal and gas clean-up [108]. The improvements also focused on identifying site specificity and a preferable well spacing. However, a few critical issues are still remaining, such as subsidence, ground water depletion and contamination, controlling ground water influx, and gas leakage. These issues associated with environmental risk are considered as the greatest threats for commercialization. Thus, the main challenge is to address the professional management of the environmental risks [108,109].
Several field trials and a few commercial applications of UCG indicate that UCG is very site -specific [10,11]. The technical and commercial risk factors can be decreased considering appropriate site selection by assessing the site, particularly the depth of the coal seam as well as of the overlying and underlying strata, and the geotechnical factors such as the strength, jointing, and deformability of the overlying strata, the surrounding geology, and the hydrogeology [10,109]. In order to proceed with site selection, more detailed knowledge is required for subsidence and water influx. There is still a lack of knowledge of the factors regarding underground water influx and subsidence. However, a more detailed understanding of modeling including the factors for subsidence and water influx can predict the conditions from the output.
Beside the technical and environmental concerns, public perception is very decisive in permitting the process for UCG [109]. Thus, cultivating public acceptance of the technology is another challenge. Recently, the UCG association has been working as an information source of UCG and they are promoting the commercial, social, and environmental benefits of UCG, providing a public and independent information service on aspects of UCG and representing it at the highest levels of governments, regulatory bodies, licensing authorities, energy groups, decision-makers, environmental groups, and the media, maintaining strong links with academia and focusing on education, training, and regulation [110]. As a result, stakeholders are gaining confidence in the new technology and the worldwide professional as well as non-professional communities are getting well educated about UCG technology. They are recognizing UCG as the clean energy source of the future with its outstanding advantages, such as the elimination of coal stock piles and coal transport and much of the disturbance at the surface, low dust and noise levels, the absence of health and safety concerns relating to underground workers, the avoidance of ash handling at power stations, and the elimination of SO2 and NOx emissions [111].
Finally, the economics of production and utilisation of the gas plays an important role for commercialization. The primary use of syngas obtained from UCG is to generate electricity; however, in order to achieve an expanded and diversified economy, the syngas are now being used as in many applications, such as the production of chemicals (e.g., hydrogen, ammonia, methanol, ethanol, urea, naptha) or liquid fuels (gasoline, diesel, jet fuel, LPG, etc.). In addition, research shows that UCG- IGCC (integrated gasification combined cycle) power plants are cheaper and carbon-efficient [112]. According to the analysis of Brand [113], the cost of syngas production is highly competitive against the available substitute on the market of Central East Europe. The scenario is different for other countries where natural gas is abundant and cheap. The most vital economic factor that delays the commercialization of UCG is the market force—lower natural gas prices. Thus, it is expected that the development of UCG will be in those countries where there is a scarcity of natural gas and an abundance of coal such as in South Africa, China, and Australia, where these developments are taking place with more vigour.
However, at present, the lack of investor confidence and the slower uptake of UCG in the economic climate are due to the absence of the information of a fully operational plant using up-to-date UCG technology [109]. As the worldwide pilot plants of some renowned companies including Linc Energy Ltd., Cougar Energy Ltd. and Carbon Energy Ltd. are working toward the commercialization of UCG and are very close to the commercialization [27], the situation is expected to change very soon in the near future.

4. UCG-CCS (carbon capture and storage)

Currently, sequestration of CO2 in a subsurface is becoming popular worldwide because of moving towards a low-carbon future. The unminable coal seams are currently being considered as carbon dioxide sinks [111,114]. The accessibility of carbon dioxide to the pores of unworked deep coal seams is limited by the porosity and permeability of coals that tend to decrease with increasing overburden pressure [115]. However, the void created due to gasification and the amount of coal chars left behind after the UCG process triggers the possibility of storing carbon dioxide in that depleted space known as UCG-CCS [116].
Considering the chemistry of gas storage, studies show that gas storage in coal occurs in both the adsorbed phase and the free gas phase. At low pressure, the carbon dioxide storage potential in coal seams is mainly a function of its adsorption capacity; however, at high pressure, it is a strong function of the accessibility to CO2 with respect to pore space for free gas [115,116,117]. Although the volume of pores in coal is much smaller than in other reservoir rocks, the storage potential of CO2 in coal significantly differs from other rocks and exceeds its open pore volume by an order of magnitude predominantly due to the surface adsorption of micropores [117]. Saghafi et al. [117] reported that coal micropore surface area can reach several hundreds of square meters per gram of solid, which makes large areas available for gas adsorption. However, the development of micropores is a function of coal rank, temperature, and coal composition [117].
There are a number of reports concerning CO2 storage potential in underground coal seams of different ranks [118]; however, there were inadequate reports of CO2 storage potential in the coal seams subsequent to gasification. Recently, case studies of the feasibility for UCG-CCS storage have been carried out in the Powder River Basin of Wyoming, US [7,93,113,119], at a selected site in Bulgaria [91,93], and at the Williston Basin, North Dakota, US [108]. Roddy and Younger [111] studied the potential of carbon storage and addressed surface subsidence and aquifer contamination as the main challenges for UCG-CCS application; however, in a later paper, Younger [120] studied hydrogeological and geomechanical aspects of UCG and proposed an action plan in order to avoid subsidence and aquifer contamination in the life cycle of UCG-CCS operations. In order to predict these risks, geomechanical and hydrogeological models are required. The models that dealt with surface subsidence and water hydrology in the UCG process can be found in recent references: Yang et al. [97], Sheng et al. [98]. From the studies, it appears that careful planning, appropriate site selection, adequate engineering and modeling, and professional management are the keys in order to implement UCG-CCS schemes successfully. It was also understood from the literature that the cavity should be located deeper than 800 m, so that CO2 can be stored in the supercritical state, allowing significantly higher utilization of the pore space [7,111]. Another interesting feature of UCG-CCS is the cavity temperature during CO2 injection. As the temperature in the cavity of the UCG process reaches more than 1000 °C, the reactor cool-down process is essential for the preparation of UCG cavities for the subsequent CO2 storage. The cooling of the cavity can be achieved either naturally very slowly; however, forced cooling with water flushing can significantly decrease the duration of the cool-down process by a factor of more than three hundred [99]. Finally, the implementation of UCG-CCS must be dependent on economic feasibility. Although very few studies are available on the economic assessment for the evaluation of UCG-CCS schemes [98,121], much remains to be done.

5. Conclusions

A number of significant modeling efforts and the outcomes of these models are presented in this review. The effect of operating conditions, such as temperature, pressure, injected gas composition and flow rate, and steam/O2 ratio on gas composition or cavity growth rate are revealed by a number of models. For most of the cases, the results are in qualitative agreement among the models. However, it is hard to compare the results quantitatively when the approach of the models and assumptions are different.
The packed bed models are applicable for highly permeable porous media only. Due to the assumption of highly porous media, Darcy’s law has been considered for fluid flow. Thus, the possibility of turbulent flow at any location of the gasification channel was ignored. None of the models considered radiation heat transfer and heat loss to the surroundings. However, most of the packed bed models considered a detailed physical and chemical model, based on experimental and theoretical comparisons and correlations. Comprehensive drying and pyrolysis sub-models were also considered in most of the packed bed models. As a result, the packed bed models exhibit good agreement for gas compositions measured from laboratory experiments. In order to calculate the heat recovery and gas composition, the model can be very effective. However, the cavity shape could not be obtained by any of the existing packed bed models due to their 1D framework. Another difficulty is to incorporate thermo-mechanical failure into the models due to their approach. As a result, the extension of these models to field trials is formidable.
The channel models overcome the limitation of the packed bed models in regards to calculating the cavity shape and size. In this method, a permeable channel between the injection and production holes is preferred. Coal is considered to be gasified only at the perimeter of the expanding permeable channel. Natural convection is considered the most important phenomenon for the mixing of the injected blast gas and the gas coming from the channel wall. All models of heat transfer were considered in most of the channel models. Some models successfully incorporated water influx. Fluid flow was not limited to Darcy’s law. As a result, in some models, the transport of heat and species inside the cavity was considered to be governed by empirical correlations for turbulent transport phenomena in enclosures. However, most of the channel models neglected drying and pyrolysis, which are considered very important in other approaches of UCG models. Very few models included thermo-mechanical failure. Despite this fact, most of the channel models were validated with the data from field trials. The channel model is found to better calculate sweep efficiency. Most of the existing channel models are in 2D. However, recently, a few 3D models have been developed from which the cavity shape and size are visibly obtained.
The coal slab models describe the process by movement of various defined regions in a coal slab, usually perpendicular to the flow of the injected blast gas. The special feature of the coal slab model is the tracking of drying and combustion front movement. These models can successfully demonstrate the drying and devolatilization behaviour of large coal particles; however, these models are yet to be validated using UCG field trial data. Due to the assumption of semi-infinite coal slab, it can be speculated that this model is only good for thick coal layers. All the models so far were developed by considering that coal slabs are one-dimensional (1D). As a result, information on cavity shape cannot be obtained from these models.
Regardless of the approach of the models, there is not a single model that includes all the important physical and chemical models for the successful prediction of UCG. Although thermo-mechanical failure is a fact for field trials, there are only very few models that consider this phenomenon with unrealistic and simple assumptions. Because of computational time, most of the earlier models are either in a 1D or 2D framework. Even the 3D models developed thus far are very simple and did not consider the details of the physical and chemical models. Considering all these facts about the UCG models, one should follow the approaches of the models depending upon one’s objective and available data at hand; for example, one may attempt to follow the approach of the channel model if he/she is interested in modeling the cavity growth, where others may be interested in the packed bed model for gas composition using fairly complicated chemical process and design parameters.
Developing a successful UCG model requires a sound understanding of physical theory as well as knowledge of details to interpret field test and laboratory experimental results in order to successfully predict the field test data and predict the performance of a new site before proceeding with gasification. It is established that each potential UCG site will be unique, so design parameters (i.e., coal seam depth, thickness, inclination, coal properties, moisture content, water influx, etc.) are needed to be dealt with first, before operating and performance parameters. That is why a number of experimental laboratory data are needed in order to develop the viable correlation of various properties to use in the model to enhance the accuracy and range of the applicability of the model. A detailed analysis is required for a specific UCG site for seam geology, properties, overburden properties, and hydrology in order to anticipate the possible types of changes in the state of the roof rock so that the approach of modeling can be carefully chosen. The future focus should be the development of a model in a three-dimensional frameworks including all the factors that affect the performance of the UCG.
Commercialization of UCG technology is a vital issue due to some challenges, primarily the environmental concerns and the low energy prices. However, due to recent repeated trial success, progress in handling the environmental issues and the improvement in costs and, finally, the global need for new clean energy sources, UCG is becoming very appealing to the energy sectors. As a result, a number of worldwide renowned companies are working effectively towards commercialization. Considering all these facts, it can be said that UCG is now on the platform close to commercialization. In addition, recent feasibility studies on UCG-CCS schemes added a low-carbon future to UCG.
A number of modeling approaches to gasification and their applicability and relevance have been carefully reviewed. Research in this area is necessary to proceed in order to able to successfully develop and implement the UCG process.

Acknowledgments

We would like to thank the Canadian Centre for Clean Coal/Carbon and Mineral Processing Technologies (C5MPT), University of Alberta, Canada for the financial support of this work.

Author Contributions

M.M. Khan conducted the study and wrote the paper that was started by A.S. Shirazi. Mmbaga also contributed broadly in defining the framework of this paper and revising the manuscript. All other authors contributed by providing their thoughts and suggestion. Finally, the work was carried out under active and direct supervision of Gupta.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

Unit
a i j ,   a k j Stoichiometric matrix
CMolar concentration of gas speciesmol/(m2s)
C i Molar concentration of species i (gas)mol/(m2s)
CwMolar concentration of watermol/(m2s)
C 2 Inertial resistance factor1/m
C p g Specific heat of gasJ/(mol K)
C p s Specific heat of solidJ/(mol K)
C p w Specific heat of waterJ/(mol K)
C 1 ε ,   C 2 ε Constant for turbulence model
D Solid thermal ConductivityW/(m K)
D i Diffusion coefficient of gas speciesm2/s
DeffEffective diffusivitym2/s
dPipe diameterM
e Specific internal energy of flowing gasJ/mol
e s Solid specific energyJ/mol
EActivation EnergyJ/mol
GNumber of gas species
GrGrashof number
G k t ,   G b Source terms for turbulence modelkg/(m·s3)
h Specific enthalpy of flowing gasJ/mol
hiSpecific enthalpy of gas species iJ/mol
hkSpecific enthalpy of solid species kJ/mol
h T Convective heat transfer coefficientW/(m2·K)
H 1 ,   H 2 Coal seam length and height respectivelyM
HiHeat of combustion of gas species iW/m3
Δ H i Heat of reaction of species iW/m3
Δ H v Heat of evaporation of waterJ/mol
jiTotal flux of gas species imol/(m2·s)
J ¯ Gas diffusion flux vectorkg/(m2·s)
KNumber of solid/condensed species
kRate constant
koFrequency factor
k t Turbulent kinetic energyJ/m3s
K d Thermal conductivity of dried coalW/K
K w Thermal conductivity of wet coalW/K
k e f f Effective thermal conductivity calculated from thermal conductivity of solid and gasW/K
kyMass transfer coefficientmol/(m3·s)
LeEntrance lengthM
m ˙ Mass flow ratekg/s
MMoisture content of coal
MkMolecular weight of solid species kkg/mol
MwMolecular weight of waterkg/mol
p PressurePa
p Pressure in blast region outside of cavity boundary layerPa
PTotal pressurePa
pwVapor pressure of waterPa
PrParameter
Pr°Base value of the parameter
ΔprPerturbation of the value of parameter
Q g Energy release in gas due to chemical reactionsJ/m3s
Q s Solid energy source rateJ/m3s
r Radius in cylindrical coordinatesM
R j Jth reaction ratemol/(m3·s)
R c i j Rate of reaction of carbon with O2, CO2, H2O and H2mol/(m3·s)
RTTotal chemical rate of reactionmol/(m3·s)
RcEffective chemical rate of reactionmol/(m3·s)
RmMass transfer rate associated with reactionmol/(m3·s)
RGas constantJ/(mol K)
rateiReaction ratekg/(m3·s)
ReReynolds number
S m Gas phase mass source due to chemical reactionskg/(m3·s)
S Specific surface area1/m
ScSchmid number
ShSherwood number
S(pr)Sensitivity of parameter
S s Solid species mass source per unit volume of bedkg/(m3·s)
SφSource term due to heterogeneous chemical reactionskg/(m3·s)
S φ Source term for porosity, calculated based on reactions in coal
sUpstream end velocitycm/min
sbDownstream end velocitycm/min
T Gas temperatureK
T s Solid temperatureK
T d Temperature in the dried core regionK
T w Temperature in the wet core regionK
T v a p Evaporation temperatureK
x , y Space coordinate
X   ( y , t ) Burn front locationM
x ¯ Location of drying frontM
x r z Total length of the reaction zoneM
XChar conversion
Y i mass fraction of species in gas phase
Y s i Mass fraction of solid component
YmSource term for kinetic energykg/(m·s3)
yiMole fraction of species in gas phase
ywMole fraction of water
ylLimiting reactant
uSuperficial velocity (gas)m/s
uwMobile water velocitym/s
U Velocity of uniform streamm/s
v ¯ Velocity vector (gas)m/s
| v | Velocity (scalar) of gasm/s
v s ¯ Superficial solid velocitym/s
v c Burn front velocitym/s
v c v Cavity growth ratem/s
VMVolatile matter at any moment
VM*Volatile matter content from proximate analysis
w k Mass fraction of solid species k
w c Weight fraction of carbon in charcoal
w k t ,   ˙ w ˙ ε User defined source terms
Greek Symbols
ρ Blast densitykg/m3
ρ g Density of gaskg/m3
ρ s Density of solidkg/m3
ρ k Density of kth solid specieskg/m3
ρ w Density of liquid waterkg/m3
ρ w o Uniform initial moisture densitykg/m3
φ porosity
φ l i m Upper limit of porosity
φ e Bed porosity external to particles
ΨStream function
μ g Viscosity of gasPa·s
μ t Turbulent viscosityPa·s
ε Rate of dissipation of turbulent kinetic energym2/s3
α permeabilitym2
α o Initial permeabilitym2
η Dimensionless position variable in vertical direction
ηrResistance parameter
γResistance parameter
ξ Dimensionless position variable in horizontal direction
β Weighing coefficient
σExperimental parameter in permeability-porosity equation
σ ε ,   σ k t Prandtl number for the turbulent kinetic energy and its rate of dissipation, respectively
ϕUser defined scalar
ϕvVelocity potential
θ Angle between x axis and radius in cavity
κRadiation extinction coefficient1/m
ΓDiffusion coefficient

References

  1. BP Statistical Review of World Energy—June 2007, Egham, Surrey, UK. BP Distribution Services, p. 45. Available online: http://www.bp.com/content/dam/bp-country/en_ru/documents/publications_PDF_eng/Statistical_review_2007.pdf (accessed on 20 October 2015).
  2. Couch, G.R. Underground Coal Gasification; IEA Clean Coal Centre: London, UK, Report No.: Contract No.: CCC/151; 2009; p. 129. [Google Scholar]
  3. Word Energy Council. Survey of Energy Resources, Word Energy Council: London, UK, 2007. Available online: http://ny.whlib.ac.cn/pdf/Survey_of_Energy_Resources_2007.pdf (accessed on 20 October 2015).
  4. Tsui, H.; Eu, C. Operating Concept of Circulating Fluidized Bed Gasifier from the Kinetic Point of View. Power Technol. 2003, 132, 176–183. [Google Scholar] [CrossRef]
  5. China Labour Bulletin. Deconstructing deadly details from China’s coal mine safety statistics CLB. 2006. Available online: http://www.clb.org.hk/en/content/deconstructing-deadly-details-chinas-coal-mine-safety-statistics (accessed on 15 October 2015).
  6. US Mine Safety and Health Administration. Statistics—Coal Mining Fatalities by State—Calendar Year. 2015. Available online: http://www.msha.gov/stats/charts/coalbystates.pdf (accessed on 15 October 2015). [Google Scholar]
  7. Shafirovich, E.; Varma, A. Underground Coal Gasification: A Brief Review of Current Status. Ind. Eng. Chem. Res. 2009, 48, 7865–7875. [Google Scholar] [CrossRef]
  8. Abdel-Hadi Eed, A.A.; Hsu, T.R. Computer Modeling of Fixed Bed Underground Coal Gasification Using the Permeation Method. J. Energy Resour. Tech. 1987, 109, 11–20. [Google Scholar] [CrossRef]
  9. Kariznovi, M.; Nourozieh, H.; Abedi, J.; Chen, Z. Simulation Study and Kinetic Parameter Estimation of Underground Coal Gasification in Alberta Reservoirs. Chem. Eng. Res. Des. 2013, 91, 464–476. [Google Scholar] [CrossRef]
  10. Bhutto, A.W.; Bazmi, A.A.; Zahedi, G. Underground Coal Gasification: From Fundamentals to Applications. Prog. Energy Combust. Sci. 2013, 39, 189–214. [Google Scholar] [CrossRef]
  11. Burton, E.; Friedmann, J.; Upadhye, R. Best Practices in Underground Coal Gasification; Lawrence Livermore National Laboratory: California, USA, 2006. Available online: http://www.purdue.edu/discoverypark/energy/assets/pdfs/cctr/BestPracticesinUCG-draft.pdf (accessed on 13 August 2015).
  12. Nourozieh, H.; Kariznovi, M.; Chen, Z.; Abedi, J. Simulation Study of Underground Coal Gasification in Alberta Reservoirs: Geological Structure and Process Modeling. Energy Fuels 2010, 24, 3540–3550. [Google Scholar] [CrossRef]
  13. Anthony, D.B.; Howard, J.B. Coal Devolatilization and Hydrogastification. AIChE J. 1976, 22, 625–656. [Google Scholar] [CrossRef]
  14. Yang, L. Numerical Study on the Underground Coal Gasification for Inclined Seams. AIChE J. 2005, 51, 3059–3071. [Google Scholar] [CrossRef]
  15. Perkins, G. Mathematical Modelling of Underground Coal Gasification. Ph.D. Thesis, The University of New South Wales, Sydney, Australia, 2005. [Google Scholar]
  16. Klimenko, A.Y. Early Ideas in Underground Coal Gasification and Their Evolution. Energies 2009, 2, 456–476. [Google Scholar] [CrossRef]
  17. Kumar, R.N.H.; Udayakumar, D.L.; Stojcevski, A.; Maung Than Oo, A.M. Underground Coal Gasification: An alternate, Economical, and Viable Solution for future Sustainability. Int. J. Eng. Sci. Invent. 2014, 3, 57–68. [Google Scholar]
  18. Su, F.; Nakanowataru, T.; Itakura, K.-I.; Ohga, K.; Deguchi, G. Evaluation of Structural Changes in the Coal Specimen Heating Process and UCG Model Experiments for Developing Efficient UCG Systems. Energies 2013, 6, 2386–2406. [Google Scholar] [CrossRef]
  19. Gregg, D.W.; Edgar, T.G. Underground Coal Gasification. AIChE J. 1978, 24, 753–781. [Google Scholar] [CrossRef]
  20. Skafa, P.V. Underground Gasification of Coal; Gosudarstvennoe Hauchno-Tekhnicheskoe Izdatel’stvo Literatury Po Gornomu Delu, Moscow, 1960; p. 324; English transl.: Report No. UCRL-Trans 10880; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1975. [Google Scholar]
  21. Gregg, D.W.; Hill, R.W.; Olness, D.U. An Overview of the Soviet Effort in Underground Coal Gasification; Report No. UCRL-52004; Lawrence Livermore National Laboratory: Berkeley, CA, USA, 1976. [Google Scholar]
  22. Derbin, Y.; Walker, J.; Wanatowski, D.; Marshall, A. Soviet Experience of Underground Coal Gasification Focusing on Surface Subsidence. Appl. Phys. Eng. 2015, 16, 839–850. [Google Scholar] [CrossRef]
  23. Thorsness, C.B.; Britten, J.A. Underground Coal Gasification Project: Final Report; Report No. UCID-21853; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1989. [Google Scholar]
  24. Blinderman, M.S.; Saulov, D.N.; Klimenko, A.Y. Forward and Reverse Combustion Linking in Underground Coal Gasification. Energy 2008, 33, 446–454. [Google Scholar] [CrossRef]
  25. Cui, Y.; Liang, J.; Wang, Z.; Zhang, X.; Fan, C.; Liang, D.; Wang, X. Forward and Reverse Combustion Gasification of Coal with Production of High-Quality Syngas in a Simulated Pilot System or In Situ Gasification. Appl. Energy 2014, 131, 9–19. [Google Scholar] [CrossRef]
  26. Bedi, P.; Binshtok, A.; Blinderman, A.; Blinderman, M.; Elliott, G. Application of the Exergy UCG Technology in CTX Projects. In Proceedings of the 2nd Underground Coal Gasification Network Workshop, Banff, AB, Canada, 22–23 August 2012.
  27. Shirazi, A.A. CFD Simulation of Underground Coal Gasification. Master’s Thesis, University of Alberta, Alberta, Canada, 2012. [Google Scholar]
  28. Davis, B.; Mallett, C. Bloodwood Ck UCG Developments. In Proceedings of the 2nd Workshop-Underground Coal Gasification Network, Banff, AB, Canada, 22–23 August 2012.
  29. Beath, A.; Craig, S.; Littleboy, A.; Mark, R.; Mallett, C. Underground Coal Gasification: Evaluating Environmental Barriers. CSIRO Exploration and Mining Report P2004/5, Kenmore, Queensland, Australia. 2004. [Google Scholar]
  30. Shanon, M.J.; Thorsness, C.B.; Hill, R.W. Early Cavity Growth during forward Burn; Report No. UCRL-84584; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1980. [Google Scholar]
  31. Hill, R.W. Burn Cavity Growth during the Hoe Creek No. 3 Underground Coal Gasification Experiment; Report No. UCRL-85173-R-1, U.S. DOE; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1981. [Google Scholar]
  32. Yeary, D.L.; Riggs, J.B. Study of Small-scale Cavity Growth Mechanisms for Underground Coal Gasification. In Situ 1987, 11, 305–327. [Google Scholar]
  33. Poon, S.S. The Combustion Rates of Texas Lignite Cores. Masters Thesis, University of Texas, Austin, TX, USA, 1985. [Google Scholar]
  34. Daggupati, S.; Mandapati, R.N.; Mahajani, S.M.; Ganesh, A.; Mathur, D.K.; Sharma, R.K.; Aghalayam, P. Laboratory Studies on Combustion Cavity Growth in Lignite Coal Blocks in The Context of Underground Coal Gasification. Energy 2010, 35, 2374–2386. [Google Scholar] [CrossRef]
  35. Yang, L. Study on the Method of Two-Phase Underground Coal Gasification with Unfixed Pumping Points. Energy Sources 2003, 25, 917–930. [Google Scholar] [CrossRef]
  36. Prabu, V.; Jaynati, S. Simulation of Cavity Formation in Underground Coal Gasification using Borehole Combustion. Energy 2011, 36, 5854–5864. [Google Scholar] [CrossRef]
  37. Stańczyk, K.; Howaniec, N.; Smoliński, A.; Świądrowski, J.; Kapusta, K.; Wiatowski, M.; Grabowski, J.; Rogut, J. Gasification of Lignite and Hard Coal with Air and Oxygen Enriched Air in A Pilot scale Ex Situ Reactor for Underground Gasification. Fuel 2011, 90, 1953–1962. [Google Scholar] [CrossRef]
  38. Thorsness, C.B.; Hill, R.W. Coal Block Gasification: Laboratory Results and Field Plans; Report No. UCRL-85906; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1981. [Google Scholar]
  39. Upadhye, R.; Burton, E.; Friedmann, J. Science and Technology Gaps in Underground Coal Gasification; Report No. UCRL-TR-222523; Lawrence Livermore National Laboratory: Livermore, CA, USA, 2006. [Google Scholar]
  40. Gunn, R.D.; Krantz, W.B. Underground Coal Gasification: Development of Theory, Laboratory Experimentation, Interpretation, & Correlation with the Hanna Field Tests: Final report; U.S. Department of Energy: Morgantown, WV, USA, 1987. [Google Scholar]
  41. Higgins, G.H. A New Concept for in-Situ Coal Gasification; Report No. UCRL-51217; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1976. [Google Scholar]
  42. Thorsness, C.B.; Rosza, R.B. In-Situ Coal Gasification Program: Model Calculations and Laboratory Experiments; Report No. UCRL-78302; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1976. [Google Scholar]
  43. Uppal, A.A.; Bhatti, A.I.; Aamir, E.; Samar, R.; Khan, S.A. Control Oriented Modeling and Optimization of One Dimensional Packed Bed Model of Underground Coal Gasification. J. Process Control 2014, 24, 269–277. [Google Scholar] [CrossRef]
  44. Dinsmoor, B.; Galland, J.M.; Edgar, T.F. The Modeling of Cavity Formation during Underground Coal Gasification. J. Petroleum Technol. 1978, 30, 695–704. [Google Scholar] [CrossRef]
  45. Gunn, R.D.; Whitman, D.L. An In-Situ Coal Gasification Model (Forward Mode) for Feasibility Studies and Design; Report No. LERC/RI-76/2; Laramie Energy Research Center: Laramie, WY, USA, 1976. [Google Scholar]
  46. Khadse, A.N.; Qayyumi, M.; Mahajani, S.; Aghalayam, P. Model for the Underground Coal Gasification (UCG) Channel Reactor Model for the Underground Coal Gasification (UCG) Channel. Int. J. Chem. React. Eng. 2006, 4, 1–25. [Google Scholar]
  47. Winslow, A.M. Numerical Model of Coal Gasification in A Packed Bed; Report No. UCRL-77627; Lawrence Livermore National Laboratory (LLNL): Berkeley, CA, USA, 1976. [Google Scholar]
  48. Thorsness, C.B.; Grens, E.A.; Sherwood, A.A. One Dimensional Model for in Situ Coal Gasification; Report No. UCRL-52523; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1978. [Google Scholar]
  49. Campbell, J.H. Pyrolysis of Subbituminous Coal in Relation to In-Situ Coal Gasification. Fuel 1978, 57, 217–224. [Google Scholar] [CrossRef]
  50. Shirazi, A.S.; Karimipour, S.; Gupta, R. Numerical Simulation and Evaluation of Cavity Growth in In Situ Coal Gasification. Ind. Eng. Chem. Res. 2013, 52, 11712–11722. [Google Scholar] [CrossRef]
  51. Thorsness, C.B.; Kang, S.W. Further Development of a General-Purpose Packed Bed Mode for Analysis of UCG Process; Report No. UCRL-92489; Lawrence Livermore National Laboratory: Berkeley, CA, USA, 1985. [Google Scholar]
  52. Harker, J.H.; Backhurst, J.R. Fuel and Energy (Energy Science and Engineering); Academic Press: London, UK, 1981. [Google Scholar]
  53. Gadsby, J.; Hinshelwood, C.N.; Sykes, K.W. The Kinetics of the Reactions of the Steam-Carbon System. Proc. R. Soc. (Lond.) Series A (Math. Phys. Sci.) 1946, 187, 129–151. [Google Scholar] [CrossRef]
  54. Lewis, W.K.; Gilliland, E.R.; Paxton, R.R. Low-Temperature Oxidation of Carbon. Ind. Eng. Chem. 1954, 46, 1327–1331. [Google Scholar] [CrossRef]
  55. Batenburg, W.D.; Biezen, N.J.; Bruining, J. A New Channel Model for Underground Coal Gasification of Thin, Deep Coal Seams. In Situ 1994, 18, 419–451. [Google Scholar]
  56. Kuyper, R.A.; van der Meer, T.H.; Bruining, J. Simulation of Underground Gasification of Thin Coal Seams. In Situ 1996, 20, 311–346. [Google Scholar]
  57. Eddy, T.L.; Schwartz, S.H. A Side Wall Burn Model for Cavity Growth in Underground Coal Gasification. J. Energy Resour. Tech. 1983, 105, 145–156. [Google Scholar] [CrossRef]
  58. Luo, Y.; Coertzen, M.; Dumble, S. Comparison of UCG Cavity Growth with CFD model predictions. In Proceedings of the 7th International Conference on CFD in the Minerals and Process Industries CSIRO, Melbourne, Australia, 9–11 December 2009.
  59. Pirlot, P.; Pirard, J.P.; Coeme, A.; Mostade, M.A. Coupling of Chemical Processes and Flow in View of the Cavity Growth Simulation of an Underground Coal Gasifier at Great Depth. In Situ 1998, 22, 141–156. [Google Scholar]
  60. Kuyper, R.A.; van der Meer, T.H.; Hoogendoorn, C.J. Turbulent Natural Convection Flow due to Combined Buoyancy Forces during Underground Coal Gasification of Thin Seams. Chem. Eng. Sci. 1994, 49, 851–861. [Google Scholar] [CrossRef]
  61. Perkins, G.; Sahajwalla, V. Modelling of Heat and Mass Transport Phenomena and Chemical Reaction in Underground Coal Gasification. Chem. Eng. Res. Des. 2007, 85, 329–343. [Google Scholar] [CrossRef]
  62. Schwartz, S.Y.; Eddy, T.L.; Mehta, K.H.; Lutz, S.A.; Binaie-Kondolojy, M.B. Cavity growth mechanisms in UCG with side wall burn gasification. In Proceedings of the SPE Annual Fall Technical Conference and Exhibition, Houston, TX, USA, 1–3 October 1978.
  63. Daney, D.E. Turbulent Natural Convection of Liquid Deuterium Hydrogen, and Nitrogen within Enclosed Vessels. Int. J. Heat Mass Transf. 1970, 19, 431–441. [Google Scholar] [CrossRef]
  64. Bailey, A.C.; Bartlett, F.; Boswinkel, H.H.; Bruining, J.; Chandelle, V.; Del Amor, G.; Garstang, J.H.; Hewel-Bundermann, H.; Hewing, G.; Krabiell, K.; et al. The Future Development of Underground Coal Gasification in Europe; A Comprehensive Report to CEC; Commission of European Communities: Brussels, Belgium, 1989. [Google Scholar]
  65. Wilks, I.H.C. The Cavity Produced by Gasifying Thin Deep Seams. In Proceedings of the 9th Underground Coal Gasification Symposium US, Washington, DC, USA, 7 August 1983; pp. 314–322.
  66. Cena, R.J.; Britten, J.A.; Thorsness, C.B. Excavation of the Partial Seam CRIP Underground Coal Gasification Test Site; Report No. 97245; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1987. [Google Scholar]
  67. Bruining, J.; Herkstroter, N.M.; Lankhorst, A.M.; Visser, C.P.; Ronde, H. Free Convection Enhanced Mass Transfer during Underground Tunnel Gasification of Coal. In Proceedings of the Thirteenth Annual Underground Coal Gasification Symposium; U.S. Department of Energy: Washington, DC, USA, 1987; pp. 225–233. [Google Scholar]
  68. Hill, R.W.; Thorsness, C.B. Summary Report on Large Block Experiments in Underground Coal Gasification, Tono Basin, Washington: Volume 1. Experimental Description and Data Analysis; Report No. UCRL-53305; Lawrence Livermore National Laboratory (LLNL) Report: Berkeley, CA, USA, 1982. [Google Scholar]
  69. Oliver, R.L.; Mason, G.M.; Spackman, L.K. Field and Laboratory Results from the TONO I (CRIP) UCG Cavity Excavation Project, WIDCO Mine Site, Centralia, Washington. Fuel Sci. Tech. 1989, 7, 1059–1120. [Google Scholar] [CrossRef]
  70. Launder, B.E.; Sharma, B.I. Application of the Energy Dissipation Model of Turbulence to the Calculation of the Flow near a Spinning Disc. Lett. Heat Mass Transfer. 1974, 1, 131–138. [Google Scholar] [CrossRef]
  71. Biezen, E.N.J. Modeling Underground Coal Gasification. Ph.D. Thesis, Delft University of Technology, Delft, The Netherlands, 1996. [Google Scholar]
  72. Saulov, D.N.; Plumb, O.A.; Klimenko, A.Y. Flame Propagation in a Gasification Channel. Energy 2010, 35, 1264–1273. [Google Scholar] [CrossRef]
  73. Chernyshev, A.B. Izbrannye Trudy; Izdatel’stvo: Moscow, USSR, 1956. [Google Scholar]
  74. Tsang, T.H.T. Modeling of Heat and Mass Transfer during Coal Block Gasification. Ph.D. Thesis, University of Texas at Austin, Austin, TX, USA, 1980. [Google Scholar]
  75. Westmoreland, P.R.; Forrester, R.C., III. Pyrolysis of Large Coal Blocks: Implications of Heat and Mass Transport Effects for In-Situ Gasification. ACS Div. Fuel Chem. 1977, 22, 93. [Google Scholar]
  76. Perkins, G.; Sahajwalla, V. A Numerical Study of the Effects of Operating Conditions and Coal Properties on Cavity Growth in Underground Coal Gasification. Energy 2006, 20, 596–608. [Google Scholar] [CrossRef]
  77. Massaquoi, J.G.M.; Riggs, J.B. Mathematical Modeling of Combustion and Gasification of a Wet Coal Slab—I: Model Development and Verification. Chem. Eng. Sci. 1983, 38, 1747–1756. [Google Scholar] [CrossRef]
  78. Massaquoi, J.G.M.; Riggs, J.B. Mathematical Modeling of Combustion and Gasification of a Wet Coal Slab—II: Mode of Combustion, Steady State Multiplicities and Extinction. Chem. Eng. Sci. 1983, 38, 1757–1766. [Google Scholar] [CrossRef]
  79. Park, K.Y.; Edgar, T.F. Modeling of Early Cavity Growth for Underground Coal Gasification. Ind. Eng. Chem. Res. 1987, 6, 237–246. [Google Scholar] [CrossRef]
  80. Perkins, G.; Sahajwalla, V. A Mathematical Model for the Chemical Reaction of a Semi-infinite Block of Coal in Underground Coal Gasification. Energy Fuels 2005, 19, 1679–1692. [Google Scholar] [CrossRef]
  81. Forrester, R.C. Two-Dimensional Studies of Coal Pyrolyis Preliminary Results. In Proceedings of the 5th Annual Underground Coal Conversion Symposium; U.S. Department of Energy: Washington, DC, USA, 1979. [Google Scholar]
  82. Wellborn, T.A. Linear Burning Rates of Texas Lignite. Master’s Thesis, University of Texas, Austin, TX, USA, 1981. [Google Scholar]
  83. Tsang, T.H.; Edgar, T.F. Modeling of Drying and Pyrolysis during Underground Coal Gasification—I. In Situ 1983, 7, 237–264. [Google Scholar]
  84. Hsia, S.J.; Wellborn, T.A.; Edgar, T.F. Measurement of Combustion Rates of Texas Lignite with Application to Underground Coal Gasification. In Proceedings of the Fourth UCC Symposium, CO, USA, 17–20 July 1978.
  85. Poon, S.S.M.; Edgar, T.F.; Park, K.Y. The Combustion Rates of Texas Lignite Cores. In Situ 1986, 10, 109–143. [Google Scholar]
  86. Chang, H.L.; Himmelblau, D.M.; Edgar, T.F. A Sweep Efficiency Model for Underground Coal Gasification. In Situ 1985, 9, 185–221. [Google Scholar]
  87. Debelle, B.; Malmendier, M.; Mostade, M.; Pirard, J.P. Modelling of Flow at Thulin Underground Coal Gasification Experiments. Fuel 1992, 71, 95–106. [Google Scholar] [CrossRef]
  88. Thorsness, C.B. Steam Tracer Experiment at the Hoe Creek No.3 Underground Coal Gasification Field Test; Report No. UCRL-53082; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1980. [Google Scholar]
  89. Pirard, J.P.; Brasseur, A.; Coeme, A.; Mostade, M.; Pirlot, P. Results of the Tracer Tests during the El Tremedal Underground Coal Gasification at Great Depth. Fuel 2000, 79, 471–478. [Google Scholar] [CrossRef]
  90. Daggupati, S.; Mandapati, R.N.; Mahajani, S.M.; Ganesh, A.; Pal, A.K.; Sharma, R.K.; Aghalayam, P. Compartment Modeling for Flow Characterization of Underground Coal Gasification Cavity. Ind. Eng. Chem. Res. 2011, 50, 277–290. [Google Scholar] [CrossRef]
  91. Batenburg, W.D.; Bruining, J. The Efficiency of Filtration Gasification. In Situ 1993, 17, 413–437. [Google Scholar]
  92. Biezen, E.N.; Bruining, J.; Molenaar, J. An Integrated 3D Model for Underground Coal Gasification. In Proceedings of the Annual Technical Conference of the Society of Petroleum Engineers, Dallas, TX, USA, 22–25 October 1995.
  93. Biezen, E.N.; Bruining, J.; Molenaar, J. An Integrated Model for Underground Coal Gasification. In Proceedings of the Annual Technical Conference of the Society of Petroleum Engineers, New Orleans, LA, USA, 25–28 September 1994.
  94. Seifi, M.; Chen, Z.; Abedi, J. Numerical Simulation of Underground Coal Gasification using the CRIP method. Can. J. Chem. Eng. 2011, 89, 1528–1535. [Google Scholar] [CrossRef]
  95. Nitao, J.J.; Camp, D.W.; Buscheck, T.A.; White, J.A.; Burton, G.C.; Wagoner, J.L.; Chen, M. Progress on a New Integrated 3-D UCG Simulator and its Initial Application. LLNL-CONF-500551. In Proceedings of the International Pittsburgh Coal Conference, Pittsburgh, PA, USA, 12–15 September 2011.
  96. Srikantiah, S.; Samdani1, G.; Mahajani1, S.; Ganesh, A.; Aghalayam, A. 2-D Modeling of Underground Coal Gasification (UCG). In Proceedings of the Comsol Conference, Bangalore, India, 17 October 2013.
  97. Yang, D.; Sarhosis, V.; Sheng, Y. Thermal-mechanical Modelling around the Cavities of Underground Coal Gasification. J. Energy Inst. 2014, 87, 321–327. [Google Scholar] [CrossRef]
  98. Sheng, Y.; Benderev, A.; Bukolska, D.; Eshiet, K.I.; da Gama, C.D.; Gorka, T.; Green, M.; Hristov, N.; Katsimpardi, I.; Kempka, T.; et al. Interdisciplinary Studies on the Technical and Economic Feasibility of Deep Underground Coal Gasification with CO2 Storage in Bulgaria. Mitig. Adapt. Strateg. Glob. Chang. 2015, 1–33. [Google Scholar] [CrossRef]
  99. Sarhosis, V.; Yang, D.; Sheng, Y.; Kempka, T. Thermal Analysis of Underground Coal Gasification Rector Cool Down for Subsequent CO2 Storage. Energy Procedia 2013, 40, 428–436. [Google Scholar] [CrossRef]
  100. Perkins, G.; Sahajwalla, V. Steady-State Model for Estimating Gas Production from Underground Coal Gasification. Energy Fuels 2008, 22, 3902–3914. [Google Scholar] [CrossRef]
  101. Lyczkowski, R.W.; Chao, Y.T. Comparison of Stefan Model with Two-Phase Model of Coal Drying. Int. J. Heat Mass Transfer 1984, 27, 1157–1169. [Google Scholar] [CrossRef]
  102. Lyczkowski, R.W. Mechanistic Theory for Drying of Porous Media; Report No.: UCRL-52456; Lawrence Livermore National Laboratory: Livermore, CA, USA, 1978. [Google Scholar]
  103. Suuberg, E.M.; Peters, W.A.; Howard, J.B. Product Composition and Kinetics of Lignite Pyrolysis. I.E.C. Process Des. Dev. 1978, 17, 37–46. [Google Scholar] [CrossRef]
  104. Boroson, M.L.; Howard, J. Product Yields and Kinetics from Vapor Phase Cracking of Wood Pyrolysis Tars. AIChE J. 1989, 35, 120–128. [Google Scholar] [CrossRef]
  105. Juntgen, H.; van Heek, K.H. Gas Release from Coal as a Function of the Rate of Heating. Fuel 1968, 47, 103–117. [Google Scholar]
  106. Van Krevelen, D.W.; van Heerden, C.; Huntjens, F.J. Physiochemical Aspects of the Pyrolysis of Coal and Related Organic Compounds. Fuel 1951, 30, 253–259. [Google Scholar]
  107. Liu, H.; Yao, H.; Yao, K.; Chen, F.; Luo, G. Characteristics of “Three zones” during Underground Coal Gasification. Adv. Mater. Res. 2012, 524–527, 56–62. [Google Scholar] [CrossRef]
  108. Gunn, R.D. Problem Solved and Problem not Solved in UCG. Laramie Energy Research Center, Laramie, Wyoming. Preprints of Papers—American Chemical Society. Div. Fuel Chem. 1977, 22, 64–75. [Google Scholar]
  109. Green, M. Underground coal gasification: On the road to commercialisation. Mod. Power Syst. 2014, 34, 22–24. [Google Scholar]
  110. Lauder, J. Global Development of UCG. UCG Association. In Proceedings of the 2nd Underground Coal Gasification Network Workshop, Banff, AB, Canada, 22–23 August 2012.
  111. Roddy, D.J.; Paul, L.; Younger, P.L. Underground Coal Gasification with CCS: A Pathway to Decarbonising Industry. Energy Environ. Sci. 2010, 3, 400–407. [Google Scholar] [CrossRef]
  112. Blinderman, M.S.; Jones, R.M. The Chinchilla IGCC Project to Date: Underground coal gasification and Environment. In Proceedings of the Gasification Technologies Conference, San Francisco, CA, USA, 27–30 October 2002.
  113. Brand, J. Progress Feedback on Wildhorse’s Hungarian UCG Project. Wildhorse Energy. In Proceedings of the 2nd Underground Coal Gasification Network Workshop, Banff, Alberta, Canada, 22–23 August 2012.
  114. Pei, P.; Zeng, Z.; He, J. Feasibility Study of Underground Coal Gasification Combined with CO2 Capture and Sequestration in Williston Basin, North Dakota. In Proceedings of the 44th US Rock Mechanics Symposium and 5th U.S.—Canada Rock Mechanics Symposium, Salt Lake City, UT, USA, 27–30 June 2010.
  115. Kempka, T.; Fernandez-Steeger, T.; Li, D.-Y.; Schulten, M.; Schluter, R.; Krooss, B.M. Carbon Dioxide Sorption Capacities of Coal Gasification Residues. Environ. Sci. Technol. 2011, 45, 1719–1723. [Google Scholar] [CrossRef] [PubMed]
  116. Ramasamy, S.; Sripada, P.P.; Khan, M.M.; Tian, S.; Trivedi, J.; Gupta, R. Adsorption Behavior of CO2 in Coal and Coal Char. Energy Fuels 2014, 28, 5241–5251. [Google Scholar] [CrossRef]
  117. Saghafi, A.; Faiz, M.; Roberts, D. CO2 Storage and Gas Diffusivity Properties of Coals from Sydney Basin, Australia. Int. J. Coal Geol. 2007, 70, 240–254. [Google Scholar] [CrossRef]
  118. Goodman, A.L.; Campus, L.M.; Schroeder, K.T. Direct Evidence of Carbon Dioxide Sorption on Argonne Premium Coals Using Attenuated Total Reflectance-Fourier Transform Infrared Spectroscopy. Energy Fuels 2005, 19, 471–476. [Google Scholar] [CrossRef]
  119. Zamzow, K.L. Underground Coal Gasification: History, Environmental Issues, and the Proposed Project at Beluga, Alaska. Available online: http://www.groundtruthtrekking.org/Documents/UCG/UCG-KZamzow-2010.pdf (accessed on 15 October 2015).
  120. Younger, P. Hydrogeological and Geomechanical Aspects of Underground Coal Gasification and its Direct Coupling to Carbon Capture and Storage. Mine Water Environ. 2011, 30, 127–140. [Google Scholar] [CrossRef]
  121. Nakaten, N.; Schlüter, R.; Azzam, R.; Kempka, T. Development of a Techno-Economic Model for Dynamic Calculation of COE, Energy Demand and CO2 Emissions of an Integrated UCG–CCS Process. Energy 2014, 66, 779–790. [Google Scholar] [CrossRef]

Share and Cite

MDPI and ACS Style

Khan, M.M.; Mmbaga, J.P.; Shirazi, A.S.; Trivedi, J.; Liu, Q.; Gupta, R. Modelling Underground Coal Gasification—A Review. Energies 2015, 8, 12603-12668. https://doi.org/10.3390/en81112331

AMA Style

Khan MM, Mmbaga JP, Shirazi AS, Trivedi J, Liu Q, Gupta R. Modelling Underground Coal Gasification—A Review. Energies. 2015; 8(11):12603-12668. https://doi.org/10.3390/en81112331

Chicago/Turabian Style

Khan, Md M., Joseph P. Mmbaga, Ahad S. Shirazi, Japan Trivedi, Qingzia Liu, and Rajender Gupta. 2015. "Modelling Underground Coal Gasification—A Review" Energies 8, no. 11: 12603-12668. https://doi.org/10.3390/en81112331

APA Style

Khan, M. M., Mmbaga, J. P., Shirazi, A. S., Trivedi, J., Liu, Q., & Gupta, R. (2015). Modelling Underground Coal Gasification—A Review. Energies, 8(11), 12603-12668. https://doi.org/10.3390/en81112331

Article Metrics

Back to TopTop