Numerical simulation of diverse thermal in situ upgrading processes for the hydrocarbon production from kerogen in oil shale reservoirs

Abstract

We investigate the productivity and product selectivity of diverse thermal in situ upgrading processes in oil shale reservoirs. In situ upgrading processes applying the ideas of Shell In situ Conversion Process, ExxonMobil Electrofrac, and Texas A&M Steamfrac are simulated by using sector models with the assumption of symmetric patterns. In-house fully functional simulator is used, which has been developed for the kerogen pyrolysis and hydrocarbon production. In the simulation cases, sensitivity analyses to the factors having major influence on the productivity and product selectivity are conducted. The effects of the temperature of vertical heaters, the spacing of hydraulic fractures, and the position of horizontal production wells are analyzed in the applied In situ Conversion Process, Electrofrac, and Steamfrac, respectively. In the applied In situ Conversion Process cases, hydrocarbon production increases with the increasing heater temperature. In the applied Electrofrac cases, hydrocarbon production increases with the increasing spacing of hydraulic fractures, even though longer time period for the process is needed. In the applied Steamfrac cases, the case of production well located at the same depth to the injection well shows the least hydrocarbon production. Among the processes, the applied In situ Conversion Process cases show the highest weight percentage of total hydrocarbon components in the produced fluid, and the applied Electrofrac cases follow it. The applied Steamfrac cases show far lower weight percentage of hydrocarbon production than the other methods. In terms of the mass ratio of produced hydrocarbon to decomposed kerogen, the applied Steamfrac cases show the largest value among the processes by aqueous phase sweeping liquid organic phase, but they also show the huge water oil mass ratio by the continuous injection of hot water. All the applied In situ Conversion Process cases and the Electrofrac case with the short spacing of hydraulic fractures show good heating efficiency by decomposing whole kerogen in the system.

Keywords

Oil shale kerogen pyrolysis in situ upgrading compositional simulation thermal process

Introduction

Oil shale is a sedimentary rock that contains solid hydrocarbon called kerogen. Kerogen generates a number of solid and fluid components including hydrocarbons when it is heated to a high temperature above 563 K. Oil shale is a promising energy source, regarding that the rich oil shale is expected to yield 2.6 trillion barrels of oil in worldwide (Crawford et al., 2008). The success likelihood of hydrocarbon production from oil shale reservoirs is closely tied to the efficiency and the effectivity of a thermal process.

In situ upgrading processes, which can be applied to all formation depths unlike the mining and subsequent surface pyrolysis, is the practical method of heating the oil shale. Besides their applicability, in situ upgrading processes have a few advantages over the mining and subsequent surface pyrolysis, such as the easier handling of materials, less water requirements, less waste disposal, and reduced emissions of pollutants (Crawford and Killen, 2010). There still exist considerations for the application of in situ upgrading processes; excess heat should be avoided for better shale oil quality and higher recovery efficiency; natural or induced fracture systems are needed for the effective distribution of heat, because shale is not a good conductor; longer heating time is needed than the mining and subsequent surface pyrolysis; effective means for sweeping the generated hydrocarbon are needed (Crawford and Killen, 2010).

In this study, we examine the effects of factors on the productivity, product selectivity, and process efficiency of the diverse in situ upgrading technologies in the oil shale reservoir system by conducting numerical simulations. We investigate the in situ upgrading methods applying the ideas of Shell In situ Conversion Process (ICP; Vinegar, 2006), ExxonMobil Electrofrac (Symington, 2006), and Texas A&M Steamfrac (Thoram and Ehlig-Economides, 2011). From the simulations, we examine the effect of heating element geometry by comparing ICP and Electrofrac, and the effect of heating material by comparing Electrofrac and Steamfrac.

ICP method involves the heating process using multiple vertical electric heaters in a hexagonal pattern as shown in Figure 1 (Vinegar, 2006). Heat is transferred from the heaters to the oil shale formation by conduction. As kerogen decomposes, vertical production wells, which are located in the centers of the hexagons, produce converted fluids. Because heat is transferred mainly by conduction in this method, a number of heaters and producers are drilled in a small area to efficiently and timely produce the converted fluid while avoiding the cracking of liquid hydrocarbons. The pressure increases by the temperature increase and subsequent fluid expansion will be a significant driving force for the production in such systems where heaters are densely installed.

Figure 1.

Concept of ICP method.

Electrofrac method uses a hydraulic fracture propped with electrical conductive material, so that it can transfer heat to the oil shale formation linearly from the heated hydraulic fracture (Symington, 2006). Calcined petroleum coke or steel shot can be used for the electrical conductive proppant (Symington, 2006). The concept of the Electrofrac is provided in Figure 2. This method is applicable in large reservoir systems by using a multilateral horizontal well. Here, each horizontal well contains a longitudinal hydraulic fracture. The fluid generated from the kerogen decomposition is produced by the vertical wells drilled between the hydraulic fractures. Because the heat is mainly transferred into the oil shale formation by conduction as the ICP method, the dimensions and geometry of hydraulic fractures are expected to significantly affect the process efficiency and productivity.

Figure 2.

Concept of Electrofrac method.

Steamfrac method involves steam or hot water injection into the horizontal well system containing multiple vertical hydraulic fractures as shown in Figure 3 (Lee, 2014; Thoram and Ehlig-Economides, 2011). Three horizontal wells are involved in this method; one horizontal well is drilled, and hydraulic fractures are made; then two additional horizontal wells are drilled at the edges of the reservoir to produce the fluids generated from the kerogen decomposition. Injected water component presents as steam (gaseous phase) or hot liquid water (aqueous phase) according to the pressure and temperature of the reservoir system. In many cases, it presents as aqueous phase by a high reservoir pressure. Because the injected hot water can flow through the preexisting fracture network as wells as the horizontal wells and hydraulic fractures, heat is transferred by convection as well as conduction.

Figure 3.

Concept of Steamfrac method.

There are a few simulation works on the numerical simulations of ICP, Electrofrac, and Steamfrac methods. Shell has been actively conducted research to expect the heating efficiency and production behavior using the ICP method in the oil shale system of Green River Formation (Fowler and Vinegar, 2009; Shen, 2009; Vinegar, 2006). In the application of ICP method, the pattern, temperature, and spacing of heating wells have been analyzed by using a non-fractured and high permeable reservoir model (Fan et al., 2010). In the previous simulation studies of Electrofrac method, the effects of hydrostatic pressure (Symington et al., 2010), effective stress (Symington et al., 2010), thermal conductivity of proppant (Lee, 2014) have been estimated. In the previous study of Steamfrac, the effects of the number of horizontal wells and injection strategy have been analyzed (Lee, 2014). A simplified economics evaluation of each in situ upgrading method has been conducted by assuming fixed and operating cost (Lee, 2014). It has revealed that the Steamfrac method was feasible in reservoirs with large volume.

In this study, we conduct the numerical simulations by applying the ideas of ICP, Electrofrac, and Steamfrac to examine the production behavior and process efficiency, which are sensitive to the implementation details of thermal processes. In the simulations, a fully functional simulator developed to describe the production behavior and system response during the in situ upgrading processes is used (Lee et al., 2016). In the sensitivity analysis in each process, the factor affecting the productivity and heating efficiency is examined. The factors are significant to the production of phases and components by controlling chemical reactions, heat conduction, and transport properties of fluids and porous media.

Model description

Summary of the simulator

The simulator has been built on FTSim, which has the built-in capability to simulate the fluid flow in porous media (Lee et al., 2016). FTSim is a variant of the TOUGH + code of Lawrence Berkeley National Laboratory (Moridis, 2014). The all known physical and chemical phenomena of kerogen pyrolysis have been described in the developed simulator (Lee et al., 2016).

It describes the non-isothermal and multiphase–multicomponent system involving the kinetic reactions of kerogen pyrolysis. Kerogen pyrolysis is a complex process and involves a great diversity of kinetic reactions resulting in a large number of components. The simulator describes the kerogen pyrolysis model obtained from the oil shale of Green River Formation by Shell (Wellington et al., 2005). This kinetic model simplifies the kerogen pyrolysis process by six reactions as provided in Table 1. The kinetic model is composed of the reactions of kerogen decomposition, cracking of heavy oil component in gaseous phase and liquid organic phase, cracking of light oil component in gaseous phase and liquid organic phase, and coking of hydrocarbon (HC) gas. The hydrocarbon components in the reactions—heavy oil, light oil, and hydrocarbon gas—are represented by pseudo components as the mixtures of multiple hydrocarbons. In our simulation, docosane (C₂₂H₄₆), undecane (C₁₁H₂₄), and ethane (C₂H₆) are used for the heavy oil, light oil, and hydrocarbon gas, respectively. The stoichiometry of the chemical model has been adjusted by considering the molecular weights of components.

Table 1.

Formulae and parameters of the kinetic reactions of kerogen pyrolysis (Wellington et al., 2005).

Reaction	Chemical formula		Frequency factor (1/s)	Activation energy (J/mole)
Kerogen decomposition	Kerogen →	0.02691 H₂O + 0.009815 heavy oil + 0.01755 light oil + 0.0402 HC gas + 0.01049 H₂ + 0.0541 CO₂ + 0.5828 prechar	4.34 × 10⁷	161,600
Cracking of heavy oil	Heavy oil (G) → Heavy oil (O) →	0.2034 light oil + 2.1153 HC gas + 16.9188 prechar 1.8269 light oil + 0.0402 HC gas + 1.8735 prechar	9.20 × 10¹¹ 9.20 × 10¹¹	206,034 206,034
Cracking of light oil	Light oil (G) → Light oil (O) →	5.1983 HC gas 0.5125 HC gas + 11.0771 prechar	6.77 × 10¹¹ 6.77 × 10¹¹	219,328 219,328
Coking of hydrocarbon gas	HC gas (G) →	3.0 H₂ + 1.8885 char	8.87 × 10¹⁵	311,432

Including nitrogen, the simulator handles the system of 10 components and 4 phases. Ten components include seven fluid components—water, heavy oil, light oil, hydrocarbon gas, hydrogen, carbon dioxide, and nitrogen—and three solid components—kerogen, prechar, and char. Four phases include three fluid phases—aqueous, liquid organic, and gaseous phases—and one solid phase. The system is either in the single phase state or in the multiphase equilibrium state, and transition between the states occurs according to the evolution or disappearance of phases. The phase evolution and disappearance are found based on the equations of state to compute the vapor pressure of liquid phases and the solubility of gas components. In our simulation cases, solid kerogen and aqueous phase present initially, and liquid organic phase, gaseous phase, and solid products of prechar and char are generated as kerogen decomposes.

The thermophysical properties of phases and the transport properties of mobile phases are computed considering the fractions of individual components distributed in the phases. To compute the properties of individual components, diverse equations of state are used. Then the mixing rule and the effect of pressure above the vapor pressure are applied. In this process, the dynamic change of phase properties as a function of system pressure, temperature, and phase composition is accurately described. Details of the computation methods can be found from the work of Lee et al. (2016).

As kerogen decomposes, and the fluid and solid products are generated, the effective porosity and the absolute permeability of porous media change every moment. From the porosity and permeability data of the Green River Formation, the effective permeability increases from 0 to 6.908 × 10⁻¹⁵ m² as the effective porosity evolves from 0 to 0.25. The simulator computes the effective porosity by accounting for the amount of solid phase in the pores, and assigns corresponding absolute permeability. It accounts for the media deformation by applying the thermo-poro-mechanics effect, and computing the evolution of effective porosity, as the system pressure and temperature change. It also accordingly computes the increasing heat convection flow with increasing effective porosity and permeability. The simulator computes the relative permeability of fluid phases and capillary pressure in the multiphase equilibrium states by using appropriate functions for two-phase or three-phase systems (Lee et al., 2016).

The simulator computes 11 simultaneous equations of mass and energy balance per element in finite volume domains, because the system involves heat and 10 mass components. The equations of mass conservation and energy conservation are transformed into residual equations. Then the Jacobian matrix is constructed to solve the matrix equation in a fully implicit method. The matrix equation can be solved only by the fully implicit method, because of its nonlinearity due to the drastic changes of mass and energy caused by the chemical reactions. Newton–Raphson method is used to compute solutions of the matrix equations.

Many oil shale reservoirs contain preexisting fracture networks, which significantly affect the flows of heat and fluid in the system. The simulator is designed to describe the fractured porous media by using Multiple Interacting Continua (MINC) method (Pruess, 1985). It divides the primary grids into the multiple homocentric secondary grids to describe the nested fracture and matrix domains. In our study, we use the two secondary grids per primary grid and assigned fracture domain and matrix domain to the outer grid and the inner grid, respectively. Global flows of heat and flow occur in the fracture domain network, while the local flows occur between the outer fracture domain and the inner matrix domain.

Governing equations

In the simulation of in situ upgrading and subsequent hydrocarbon production, governing equations of mass conservation and energy conservation are considered as shown in Table 2. Here, chemical reactions affect the mass and heat accumulations. The mass conservation equation for each component is derived by applying the mass accumulation caused by fluid transport, sources and sinks, and mass changes due to the chemical reactions. Darcy’s law is used to compute the fluid flux term. Here,

M^{κ}

M_{react}^{κ}

, and

M_{product}^{κ}

are the mass accumulation, mass decrease of reactant, and mass increase of product, respectively.

M_{react}^{κ}

and

M_{product}^{κ}

are computed by multiplying the reaction rate to the time step considering stoichiometry. Superscript κ indicates components, and subscript β indicates phases.

F^{κ}

and

q^{κ}

indicate the flux and source/sink term, respectively.

Table 2.

Mathematical equations of the simulation.

	Governing equations	Supplementary equations
Mass conservation	$\frac{d}{dt} \int_{V_{n}} (M^{κ} + M_{react}^{κ} + M_{product}^{κ}) d V = \int_{Γ_{n}} F^{κ} \cdot n d A + \int_{V_{n}} q^{κ} d V$	$M^{κ} = \sum_{β \equiv A, O, G, S} φ S_{β} ρ_{β} X_{β}^{κ}$ $M_{react}^{κ} = \sum_{k} r_{κ, k} Δ t,$ $M_{product}^{κ} = - \sum_{k} r_{κ, k} s_{κ, k} Δ t$ $F^{κ} = \sum_{β \equiv A, O, G} F_{β} X_{β}^{κ}$ $F_{β} = - k \frac{k_{r β} ρ_{β}}{μ_{β}} (\nabla P_{β} - ρ_{β} g)$
Energy conservation	$\frac{d}{dt} \int_{V_{n}} (M^{θ} + M_{reaction}^{θ}) d V = \int_{Γ_{n}} F^{θ} \cdot n d A + \int_{V_{n}} q^{θ} d V$	$M^{θ} = (1 - ϕ) ρ_{R} \int_{T_{ref}}^{T} C_{R} d T + \sum_{β \equiv A, O, G, S} φ S_{β} ρ_{β} U_{β}$ $F^{θ} = - K \nabla T + \sum_{β \equiv A, O, G} h_{β} F_{β}$
Chemical reaction	$r_{κ, k} = K_{k} C_{β}^{κ}$	$K_{k} = A_{k} \exp (- E_{k} / (R_{g} T)),$ $C_{β}^{κ} = φ S_{β} ρ_{β} X_{β}^{κ}$

The energy conservation equation is derived by applying the heat accumulation caused by heat flow, sources and sinks of heat, and reaction heat. Here, $M^{θ}$ and $M_{reaction}^{θ}$ are the terms of heat accumulation and reaction heat. Heat accumulation term is computed by the combination of internal energies of rock and fluid phases. $F^{θ}$ and $q^{θ}$ indicate the terms of heat flux and source/sink of heat, respectively. Heat flux term is computed by adding the conductive and convective heat flows. In the heat convection term, the enthalpy of phases ( $h_{β}$ ) is considered, and the enthalpy departure is applied to account for the real gas behavior in the gaseous phase.

For the mathematical implementation of chemical reactions, the first rate law is used to compute the reaction rates as presented in the prevalent studies (Braun and Burnham, 1990, 1992). Here, $r_{κ, k}$ , K_k, and $C_{β}^{κ}$ indicate the reaction rate, reaction rate constant, and the concentration of component, respectively. The reaction rate constant is computed by using frequency factor and activation energy. The chemical reactions occur from when the system temperature approaches their activation temperatures.

Reservoir properties and initial condition

The physical, transport, and thermal properties of the reservoir in our simulations are provided in Table 3. The scaled saturation is obtained by the phase saturation divided by the total of the fluid phases (mobile phases) saturations. The initial effective porosity, oil shale grade, and fracture spacing have been obtained from the history matching with the field production data of the Green River Formation (Lee et al., 2016).

Table 3.

Reservoir properties and initial condition for the simulation cases.

Reservoir property	Value	Initial condition	Value
Rock density (kg/m³)	2000	Reservoir pressure (bar)	207
Rock thermal conductivity (W/(m·K))	2	Reservoir temperature (K)	308.15
Oil shale grade (gal/ton)	22.5	Scaled saturation of aqueous phase	1
Permeability anisotropy	0.1	Initial media porosity	0.316
Fracture spacing (m)	3	Initial effective porosity	0.04
Fracture domain permeability (m²)	1.480 × 10⁻¹³	Initial absolute permeability of matrix domain (m²)	1.228 × 10⁻¹⁶
Volume fraction of fracture domain	0.01	Initial absolute permeability of matrix domain (m²)	1.228 × 10⁻¹⁶
Irreducible scaled saturations S_irA, S_irO, S_irG* in matrix domain	0.30, 0.20, 0.05	Kerogen volumetric fraction in the bulk oil shale rock	0.276
Irreducible scaled saturations S_irA, S_irO, S_irG* in fracture domain	0.03, 0.02, 0.005

In the initial condition of the reservoir system, the fracture domain and the void space of pores were filled with aqueous phase. The kerogen volume fraction has been computed from the oil shale grade; the sum of the kerogen volume fraction and the initial effective porosity gives the initial media porosity, the total porosity including solid phase. Initial effective permeability has been obtained by inserting the initial effective porosity into the permeability model of the Green River Formation (Kibodeaux, 2014). The changing absolute permeability has been computed as the evolution of effective porosity (Nimblett and Ruppel, 2003).

Simulation models and strategies

The simulation cases applying the ideas of ICP, Electrofrac, and Steamfrac have been conducted. In each process, we assume that the principle of symmetry can be applied to the problems. A sector model, which is the minimum unit of the reservoir system, is simulated in each process by applying the no-flow of fluid and heat at the boundaries. Each model includes the overburden and underburden formations to account for the heat loss. Zero permeability is applied to the overburden and underburden formations, in order to only provide the heat loss. Each process is continued until the kerogen mass in place approached its plateau.

The simulation cases of the applied ICP have been conducted by considering the reservoir configuration as shown in Figure 4. The multiple hexagons of heaters and producers are considered; in one hexagon, 16 electric heaters and 2 producers are located. A quarter of one hexagon as a sector model is simulated as shown in Figure 5. On the assumption that the spacing and the pattern of heaters are fixed as proposed by the Shell’s MDP-S (Fowler and Vinegar, 2009), the temperature of heaters is the most important factor for the productivity and the heating efficiency. In this regard, the sensitivity analysis of production behavior to the temperature of heaters has been conducted. The temperatures of heaters = 598.15 K, 603.15 K, 608.15 K, and 613.15 K are used. It is expected that this range of temperature makes the kerogen decomposition active while avoiding the excessive cracking of hydrocarbons, especially for the cracking of liquid organic phase.

Figure 4.

Top views of the reservoir system and the quarter of one hexagon in the applied ICP cases: (a) Top view of the reservoir system and (b) Top view of the quarter of one hexagon.

Figure 5.

Simulation model of the applied ICP cases.

The top view of the reservoir system in the applied Electrofrac simulation is provided in Figure 6. A multilateral horizontal well system is considered, where the longitudinal hydraulic fracture has been made in each branch well. Each horizontal well transfers heat into its hydraulic fracture, which is propped with electrically conductive proppant. Two vertical production wells are located between the hydraulic fractures to produce converted fluids including hydrocarbons. The sector model of an eighth of one hydraulic fracture system as shown in Figure 7, which is the minimum unit of the system, is simulated. The hydraulic fractures are 30 m-high and 30 m-wide, and the temperature of each horizontal well is kept as 613.15 K during the whole heating and production processes. The sensitivity analysis of production behavior to the spacing of hydraulic fractures has been performed by using the spacing = 10 m, 15 m, and 20 m, because it was expected to significantly affect the process efficiency as well as productivity. Because the spacing of hydraulic fractures affects the size of the one fracture unit of the reservoir, the duration necessary for the heating and production processes varies with the spacing of hydraulic fractures.

Figure 6.

Top view of the reservoir system in the applied Electrofrac cases.

Figure 7.

Simulation model of the applied Electrofrac cases.

The top view of the reservoir system in the applied Steamfrac simulation is provided in Figure 8. A horizontal well system with multistage transverse vertical fractures is considered. Two production wells have been drilled at the side edges of the reservoir after fracturing. The hydraulic fractures are 30 m-high and 30 m-wide and evenly distributed with the spacing of 15 m. The minimum unit of the reservoir system is simulated as the cases of applied ICP and Electrofrac methods. The sector model of a quarter of one hydraulic fracture system is simulated as shown in Figure 9. The sensitivity analysis of production behavior to the position of the production wells has been conducted, because it was expected to significantly affect the heat distribution in the reservoir and following fluid production. In the Figure 9, the well 1 is the injection well located at the center of the reservoir, which has been used to make the hydraulic fractures, and the well 2-1, well 2-2, and well 2-3 are the production wells located at the different depths; each of the production well is used in the each case of sensitivity analysis. In the simulation cases, the injected water component presents in aqueous phase (hot water) due to the high pressure of the reservoir system. Hot water is injected into the well 1 during the whole process; injection rate per one hydraulic fracture unit is 0.4 kg/s. Simultaneously, the fluid is produced from the well 2-1, well 2-2, and well 2-3 in the cases 1, 2, and 3, respectively.

Figure 8.

Top view of the reservoir system in the applied Steamfrac cases.

Figure 9.

Simulation model and the positions of horizontal wells on the hydraulic fracture in the applied Steamfrac cases.

Results and discussion

Applied ICP cases

Kerogen mass in place and phases productions of the applied ICP cases are provided in Figure 10. Figure 10(a) shows that whole kerogen in the system has decomposed at the end of the process in every case. From the production rates of aqueous phase in Figure 10(a), it is found that the aqueous phase production rate increases from t = 45 days to t = 90, 80, 70, and 70 days, in the cases 1, 2, 3 and 4, respectively. The peaks of the aqueous phase production rates are 0.010 kg/s (=7.522 Stock Tank Barrels/Day [STB/D]), 0.010 kg/s (=7.522 STB/D), 0.012 kg/s (=8.514 STB/D), and 0.013 kg/s (=9.222 STB/D) in those cases, respectively. After the peaks, the production rates of aqueous phase gradually decrease until they approach zero, because the generation of hydrocarbons and the evaporation of aqueous phase by the increase of temperature. As the heater temperature increases, the production time and the cumulative production of aqueous phase decrease; Table 4 shows the result quantitatively. This has been caused by that more aqueous phase vaporized, and more liquid organic phase and gaseous phase were produced with increasing heater temperature.

Figure 10.

Applied ICP cases–kerogen mass in place and phases productions in one hexagon: (a) Kerogen mass in place, (b) Production rates of aqueous phase, (c) Production rates of liquid organic phase and (d) Production rates of gaseous phase.

Table 4.

Applied ICP cases–phases productions, gas oil mass ratio, and water oil mass ratio in one hexagon.

Cases	Heater temperature	Aqueous production (kg)	Liquid organic production (kg)	Gaseous production (kg)	Gas oil mass ratio (kg/kg)	Water oil mass /ratio (kg/kg)
Case 1	598.15 K	38,348 (322 STB)	67,645 (756 STB)	46,060 (1789 MSCF)	0.681	0.567
Case 2	603.15 K	37,430 (314 STB)	68,398 (765 STB)	46,581 (1809 MSCF)	0.681	0.547
Case 3	608.15 K	36,886 (309 STB)	70,157 (784 STB)	47,601 (1849 MSCF)	0.678	0.526
Case 4	613.15 K	36,341 (305 STB)	73,910 (826 STB)	48,659 (1890 MSCF)	0.658	0.492

ICP: In situ Conversion Process.

From Figure 10(c) and (d), it is found that the peaks of liquid organic phase production and gaseous phase production increase with increasing heater temperature. The peaks of liquid organic phase production rate are 7.618 × 10⁻³kg/s (=7.359 STB/D), 8.181 × 10⁻³kg/s (=7.904 STB/D), 9.491 × 10⁻³kg/s (=9.170 STB/D), and 0.011 kg/s (=10.338 STB/D); the peaks of gaseous phase production rate are 3.987 × 10⁻³kg/s (=11.060 MSCF/D), 4.332 × 10⁻³kg/s (=12.012 MSCF/D), 5.125 × 10⁻³kg/s (=14.217 MSCF/D), and 5.840 × 10⁻³kg/s (=16.199 MSCF/D) in the cases 1, 2, 3, and 4, respectively. The oscillations of production rates of liquid organic and gaseous phases have been caused by the decreasing wellbore flowing pressure of the production well and the flow of newly generated oil and gas from the decomposition of kerogen.

The quantitative cumulative productions of liquid organic phase and gaseous phase in one hexagon are found in the Table 4. The liquid organic phase production shows the highest mass among the productions of fluid phases in every case. Gas oil mass ratio—the mass ratio of gaseous phase production to the liquid organic phase production—decreases, as the heater temperature increases; this is because the rate of increase is higher in the liquid organic phase production than the gaseous phase production with increasing heater temperature. The decrease of liquid organic phase viscosity leads to more cumulative production of liquid organic phase with increasing heater temperature. Water oil mass ratio—the mass ratio of aqueous phase production to the liquid organic phase production—also decreases with increasing heater temperature, by the more production of liquid organic phase and the more evaporation of aqueous phase.

Table 5 shows the masses of decomposed kerogen and produced fluid components in one hexagon. The percentages of fluid components indicate the weight percentages of fluid components in the produced fluid. The total weight percentages of hydrocarbon components (heavy oil, light oil, and hydrocarbon gas) in the produced fluids are 65.66%, 64.73%, 63.91%, and 63.38% in the cases 1, 2, 3, and 4, respectively. The percentages of produced hydrocarbon mass to the decomposed kerogen mass are 14.50%, 13.88%, 13.41%, and 13.33% in the cases 1, 2, 3, and 4, respectively. The masses of the produced components are affected by the change of production rates of phases, which are changed by the vapor pressure of components, relative permeability of phases, and chemical reaction rates. It has caused the first increase of production of C₂H₆ and H₂ and the later decrease of production with increasing heater temperature.

Table 5.

Applied ICP cases–masses of decomposed kerogen and produced components in one hexagon.

Heater temperature	Decomposed kerogen (kg)	Produced fluid components (kg, wt%)
Heater temperature	Decomposed kerogen (kg)	C₂₂H₄₆	C₁₁H₂₄	C₂H₆	H₂	CO₂	H₂O
598.15 K	7.257 × 10⁵	3.680 × 10⁴ (24.11%)	3.758 × 10⁴ (24.62%)	2.235 × 10⁴ (14.64%)	9.744 × 10² (0.64%)	1.012 × 10⁴ (6.63%)	4.480 × 10⁴ (29.35%)
603.15 K	7.266 × 10⁵	3.692 × 10⁴ (24.22%)	3.781 × 10⁴ (24.80%)	2.271 × 10⁴ (14.90%)	9.884 × 10² (0.65%)	1.008 × 10⁴ (6.61%)	4.394 × 10⁴ (28.82%)
608.15 K	7.268 × 10⁵	3.848 × 10⁴ (24.69%)	3.955 × 10⁴ (25.38%)	2.283 × 10⁴ (14.65%)	9.942 × 10² (0.64%)	1.022 × 10⁴ (6.56%)	4.375 × 10⁴ (28.08%)
613.15 K	7.268 × 10⁵	4.072 × 10⁴ (25.37%)	4.221 × 10⁴ (26.29%)	2.247 × 10⁴ (14.00%)	9.794 × 10² (0.61%)	1.026 × 10⁴ (6.39%)	4.387 × 10⁴ (27.33%)

ICP: In situ Conversion Process.

A grid block at the center of the simulation model has been monitored to measure the kerogen volume fraction in solid phase, saturation of solid phase, effective porosity, and absolute permeability of the matrix. Figure 11 shows the observation results. In Figure 11(a), it is found that kerogen actively decomposes from 80 to 85 days in each case. In Figure 11(b), the solid phase saturation at the monitoring point gradually decreases from t = 0 to 80 days, because the saturation of aqueous phase increases due to the expansion as the reservoir pressure increases by heating. From t = 80 days, the saturation of solid phase rapidly decreases by kerogen decomposition in every case until it approaches constant value. The saturations of solid phase at the end of the simulation are 36.66%, 36.70%, 36.78%, and 37.91% in the cases 1, 2, 3, and 4, respectively. The slight increase of the solid phase saturation is caused by the increase of the masses of solid products as the heater temperature increases, by the more vigorous chemical reactions.

Figure 11.

Observation results at the monitoring point in the applied ICP cases: (a) Kerogen volume fraction, (b) Saturation of solid phase, (c) Effective porosiy and (d) Absolute permeability of matrix domain.

Figure 11(c) and (d) shows the evolution of effective porosity and absolute permeability of the matrix domain as the process of heating and production proceeds. From t = 0 to 80 days, the effective porosity and corresponding absolute permeability gradually increase by increasing reservoir pressure due to the heating, and following fluid expansion. As kerogen decomposes, the effective porosity and the absolute permeability rapidly increase after t = 80 days.

At the end of the simulation, the values of effective porosity are 20.06%, 20.05%, 20.03%, and 19.99% in the cases 1, 2, 3, and 4, respectively. The values of corresponding absolute permeability of the matrix are 3.090 × 10⁻¹⁵ m², 3.086 × 10⁻¹⁵ m², 3.079 × 10⁻¹⁵ m², and 3.066 × 10⁻¹⁵ m², respectively. The media porosity including solid phase at the end of the simulation is around 31.7% in every case, and the presence of solid phase reduces the absolute permeability of the matrix. The decrease of absolute permeability due to the presence of solid phase is 4.736 × 10⁻¹⁵ m².

Applied Electrofrac cases

Kerogen mass in place and phases productions of the applied Electrofrac cases are provided in Figure 12. Figure 12(a) shows that the kerogen mass in place stops decreasing and reaches to a certain value at the end of the process in each case. Almost whole kerogen has decomposed in the case of hydraulic fracture spacing = 10 m, while there has remained kerogen mass in place in the cases of hydraulic fracture spacing = 15 m and 20 m. Figure 12 shows the phases productions of the simulation cases. In Figure 12(b), it is found that the aqueous phase production rate reaches to a high value at the beginning of the process and gradually decreases in each case. This is because the distance between the hydraulic fracture and the production well is far, and only aqueous phase is produced at the early time of the process. The peaks of the aqueous phase production rates are 7.429 × 10⁻³kg/s (=5.383 STB/D), 7.043 × 10⁻³kg/s (=5.103 STB/D), and 7.095 × 10⁻³kg/s (=5.141 STB/D) in the cases 1, 2, and 3, respectively. The unstable aqueous phase production rate has been caused by the compositional effects of decreasing wellbore flowing pressure of production well, continuous decomposition of kerogen into oil, gas, water, and other components, and changes of relative permeability of phases.

Figure 12.

Applied Electrofrac cases–kerogen mass in place and phases productions in one hydraulic fracture unit: (a) Kerogen mass in place, (b) Production rates of aqueous phase, (c) Production rates of liquid organic phase and (d) Production rates of gaseous phase.

In Figure 12(c) and (d), it is found that the liquid organic phase production and gaseous phase production peak at the later time than the aqueous phase production. The peaks of liquid organic phase production rate are 5.989 × 10⁻³kg/s (=5.786 STB/D), 3.317 × 10⁻³kg/s (=3.205 STB/D), and 2.537 × 10⁻³kg/s (=2.451 STB/D); the peaks of gaseous phase production rate are 8.505 × 10⁻³kg/s (=23.590 Thousand standard cubic feet per day [MSCF/D]), 4.624 × 10⁻³kg/s (=12.827 MSCF/D), and 3.822 × 10⁻³kg/s (=10.602 MSCF/D) in the cases 1, 2, and 3, respectively. As the spacing of hydraulic fractures increases, the peaks of the production rates of phases become smaller, but the more amount of phases is produced by the longer time period of the process and the larger reservoir size. The quantitative amount of the cumulative productions is provided in Table 6. The aqueous phase production shows the highest mass among the produced fluid phases in every case.

Table 6.

Applied Electrofrac cases–phases productions, gas oil mass ratio, and water oil mass ratio in one fracture unit.

Cases	Spacing of hydraulic fractures	Aqueous production (kg)	Liquid organic production (kg)	Gaseous production (kg)	Gas oil mass ratio (kg/kg)	Water oil mass ratio (kg/kg)
Case 1	10 m	1.992 × 10⁵ (1671 STB)	1.204 × 10⁵ (1346 STB)	2.052 × 10⁵ (7,976 MSCF)	1.704	1.654
Case 2	15 m	3.096 × 10⁵ (2597 STB)	1.451 × 10⁵ (1622 STB)	2.362 × 10⁵ (9174 MSCF)	1.628	2.134
Case 3	20 m	4.155 × 10⁵ (3485 STB)	1.695 × 10⁵ (1895 STB)	2.923 × 10⁵ (11,354 MSCF)	1.724	2.451

Table 7 shows the masses of decomposed kerogen and produced fluid components in one fracture unit at the end of the simulations. From the weight percentages of produced fluid components, it is found that the water component took up the highest percentage, and the light oil, hydrocarbon gas, and heavy oil components follow it. Comparing with the applied ICP cases, lower weight percentage of oil components is produced. This is because that the formation evenly approaches the high temperature in the ICP cases at the early time by the many number of wells in the small area of the reservoir, while the formation temperature reaches to the high temperature slowly in the applied Electrofrac cases. It leads to the higher viscosity and lower mobility of liquid organic phase in these cases.

Table 7.

Applied Electrofrac cases–masses of decomposed kerogen and produced components in one fracture unit.

Spacing of hydraulic fractures	Decomposed kerogen (kg)	Produced fluid components (kg, wt%)
Spacing of hydraulic fractures	Decomposed kerogen (kg)	C₂₂H₄₆	C₁₁H₂₄	C₂H₆	H₂	CO₂	H₂O
10 m	2.536 × 10⁶	6.232 × 10⁴ (11.86%)	1.013 × 10⁵ (19.28%)	7.417 × 10⁴ (14.12%)	3.232 × 10³ (0.62%)	3.535 × 10⁴ (6.73%)	2.490 × 10⁵ (47.39%)
15 m	3.254 × 10⁶	7.109 × 10⁴ (10.29%)	1.147 × 10⁵ (16.60%)	9.783 × 10⁴ (14.16%)	4.262 × 10³ (0.62%)	4.636 × 10⁴ (6.71%)	3.567 × 10⁵ (51.62%)
20 m	4.164 × 10⁶	8.113 × 10⁴ (9.25%)	1.361 × 10⁵ (15.51%)	1.265 × 10⁵ (14.42%)	5.514 × 10³ (0.63%)	6.007 × 10⁴ (6.85%)	4.680 × 10⁵ (53.34%)

The total weight percentages of hydrocarbon components in the produced fluids are 42.56%, 41.05%, and 39.18% in the cases 1, 2, and 3, respectively. The percentages of produced hydrocarbon mass to the decomposed kerogen mass are 9.38%, 8.72%, and 8.26% in the cases 1, 2, and 3, respectively.

Applied Steamfrac cases

Kerogen mass in place and phases productions of the applied Steamfrac cases are provided in Figure 13. Figure 13(a) shows that the more kerogen has remained in the system after the process in the case 2, where the production well is located at the same depth to the injection well. In Figure 13 and Table 8, it is found that the total hydrocarbon production in the case 2 is smaller than that of the case 1 or case 3. The kerogen distribution profiles of each case in Figure 14 provide the reason; the volumes of kerogen decomposition are larger in the case 1 and case 3 than in the case 2 at the end of the process.

Figure 13.

Applied Steamfrac cases–kerogen mass in place and phases productions in one hydraulic fracture unit: production well depth = 7.5 m above the injection well (case 1), production well depth = same depth to the injection well (middle of the reservoir) (case 2), production well depth = 7.5 m below the injection well (case 3): (a) Kerogen mass in place, (b) Production rates of aqueous phase, (c) Production rates of liquid organic phase and (d) Production rates of gaseous phase.

Table 8.

Applied Steamfrac cases–phases productions, gas oil mass ratio, and water oil mass ratio in one fracture unit: production well depth = 7.5 m above the injection well (case 1), production well depth = same depth to the injection well (middle of the reservoir) (case 2), production well depth = 7.5 m below the injection well (case 3).

Cases	Aqueous production (kg)	Liquid organic production (kg)	Gaseous production (kg)	Gas oil mass ratio (kg/kg)	Water oil mass ratio (kg/kg)
Case 1	3.417 × 10⁷ (2.866 × 10⁵ STB)	5.602 × 10⁵ (6264 STB)	2.251 × 10⁴ (874 MSCF)	0.040	60.996
Case 2	3.446 × 10⁷ (2.890 × 10⁵ STB)	3.452 × 10⁵ (3860 STB)	5.003 (194 MSCF)	0.014	99.826
Case 3	3.419 × 10⁷ (2.868 × 10⁵ STB)	5.785 × 10⁵ (6469 STB)	1.939 × 10⁴ (753 MSCF)	0.034	59.101

Figure 14.

Applied Steamfrac cases–profiles of kerogen volume fraction in solid phase at t = 1000 days: production well depth = 7.5 m above the injection well (case 1), production well depth = same depth to the injection well (middle of the reservoir) (case 2), production well depth = 7.5 m below the injection well (case 3): (a) X = 0.04 m (case 1), (b) X = 3.75 m (case 1), (c) X = 7.5 m (case 1), (d) X = 0.04 m (case 2), (e) X = 3.75 m (case 2), (f) X = 7.5 m (case 2), (g) X = 0.04 m (case 3), (h) X = 3.75 m (case 3) and (i) X = 7.5 m (case 3).

In Figure 13(b), it is found that the aqueous phase production rate remains the same to the injection rate of hot water in every case. The cumulative productions of liquid organic phase are higher in the case 1 and case 3 than in the case 2. The peaks of liquid organic phase production rate are 0.014 kg/s (=13.098 STB/D), 7.437 × 10⁻³kg/s (=7.185 STB/D), 0.015 kg/s (=14.063 STB/D) in the cases 1, 2, and 3, respectively. Figure 13(e) shows that the case 2 has the later start and the lower rate of gaseous phase production than the case 1 or case 3. This is caused by that the small amount of generated gaseous phase rises rather than flows to the production wells in the case 2. The case 1 has the highest cumulative production of gaseous phase and the highest gas oil mass ratio as provided in Figure 13(d) and Table 8. Comparing with the applied ICP and Electrofrac cases, this process shows the much lower gas oil mass ratio because of dominant productions of aqueous phase that sweeps liquid organic phase. The peaks of gaseous phase production rate are 8.846 × 10⁻⁴kg/s (=2.454 MSCF/D), 2.823 × 10⁻⁴kg/s (=0.783 MSCF/D), 1.545 × 10⁻³kg/s (=4.286 MSCF/D) in the cases 1, 2, and 3, respectively. As expected, the water oil mass ratio is far higher than the applied ICP and Electrofrac cases in this process.

Table 9 shows the masses of decomposed kerogen and produced fluid components in one fracture unit at the end of the process in each case. The case 3 and case 2 have the highest and the lowest mass of decomposed kerogen, respectively. From the weight percentages of produced fluid components, it is found that the water component takes up the highest percentage, and the heavy oil, light oil, and hydrocarbon gas components follow it. Comparing with the applied cases of ICP and Electrofrac, the production of water component is remarkably high in the Steamfrac cases due to the continuous hot water injection. Especially, the case 2 shows the highest weight percentage of water component production, and this turns out to be less efficient than the cases 1 and 3. The weight percentage of the produced hydrocarbon components (heavy oil, light oil, and hydrocarbon gas) to the mass of decomposed kerogen in the applied Steamfrac method is higher than that of the applied ICP and Electrofrac methods. This is because that the injected hot water provides the driving force for the converted liquid hydrocarbons to flow to the production well.

Table 9.

Applied Steamfrac cases–masses of decomposed kerogen and produced components in one fracture unit: production well depth = 7.5 m above the injection well (case 1), production well depth = same depth to the injection well (middle of the reservoir) (case 2), production well depth = 7.5 m below the injection well (case 3).

Position of production wells	Decomposed kerogen (kg)	Produced fluid components (kg, wt%)
Position of production wells	Decomposed kerogen (kg)	C₂₂H₄₆	C₁₁H₂₄	C₂H₆	H₂	CO₂	H₂O
Case 1	1.569 × 10⁶	2.911 × 10⁵ (0.84%)	2.712 × 10⁵ (0.78%)	9.090 × 10⁴ (0.26%)	2.342 × 10³ (0.007%)	2.071 × 10⁴ (0.06%)	3.408 × 10⁷ (98.05%)
Case 2	9.819 × 10⁵	1.798 × 10⁵ (0.52%)	1.662 × 10⁵ (0.48%)	3.624 × 10⁴ (0.10%)	1.392 × 10³ (0.004%)	1.440 × 10⁴ (0.04%)	3.441 × 10⁷ (98.86%)
Case 3	1.610 × 10⁶	2.990 × 10⁵ (0.86%)	2.816 × 10⁵ (0.81%)	8.438 × 10⁴ (0.24%)	2.160 × 10³ (0.006%)	2.122 × 10⁴ (0.06%)	3.410 × 10⁷ (98.02%)

The total weight percentages of hydrocarbon components (heavy oil, light oil, and hydrocarbon gas) in the produced fluids are 1.88%, 1.10%, and 1.91% in the cases 1, 2, and 3, respectively. The percentages of produced hydrocarbon mass to the decomposed kerogen mass are 41.64%, 38.93%, and 41.29% in the cases 1, 2, and 3, respectively.

Conclusion

In this study, we investigated the diverse in situ upgrading processes applying the ideas of ICP, Electrofrac, and Steamfrac methods. The effects of each process on the productions of phases and components have been analyzed. From the simulations, we have the following specific conclusions.

Applied ICP cases have showed increasing hydrocarbon production with increasing heater temperature.

Applied Electrofrac cases have showed increasing hydrocarbon production with increasing spacing of hydraulic fractures.

Applied Steamfrac cases have showed that the case of the production well at the middle of the reservoir (case 2) had the lowest hydrocarbon production due to the least distribution of heat among the cases.

In every case of applied ICP and one case of applied Electrofrac, whole kerogen in the system has decomposed.

The applied ICP cases have showed the highest total weight percentage of hydrocarbon components in the produced fluids, and the applied Electrofrac cases have followed it.

The weight percentage of produced hydrocarbons to the decomposed kerogen was highest in the cases of applied Steamfrac, by the continuously injected hot water sweeping the liquid organic phase.

The improved results can be brought by considering following issues.

Sector models have been used in the simulations, by assuming that the system responded identically in each repeatable patterns of minimum units. However, the better description of subsequent heating and boundary conditions will achieve more realistic simulations.

Coarse grids systems have been used in the simulations for the computation efficiency. Finer grids systems will reduce the numerical dispersion.

By the same reason to (2), the hydraulic fracture systems of small sizes have been used in the simulations of applied Electrofrac and Steamfrac.

Differing application strategies such as pattern geometry and heating period in each in situ upgrading process are expected to produce very variable results of productivity, product selectivity, and heating efficiency.

The simulations have revealed technical feasibility of each process. Envisioned works include the evaluation of energy efficiency and process economics of each process.

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors appreciate Crisman Institute for Petroleum Research of Texas A&M University for supporting this study.

References

Braun

Burnham

(1990) Mathematical model of oil generation, degradation, and expulsion. Energy & Fuels 4(2): 132–146.

Braun

Burnham

(1992) PMOD: A flexible model of oil and gas generation, cracking, and expulsion. Organic Geochemistry 19(1): 161–172.

Crawford P, Biglarbigi K, Dammer A, et al. (2008) Advances in world oil-shale production technologies. In: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, Denver, Colorado, 21–24 September.

Crawford

Killen

(2010) New Challenges and Directions in Oil Shale Development Technologies, ACS Symposium Series, Oxford: Oxford University Press, pp. 21–60.

Fan

Durlofsky

Tchelepi

(2010) Numerical simulation of the in-situ upgrading of oil shale. SPE Journal 15(02): 368–381.

Fowler TD and Vinegar HJ (2009) Oil shale ICP-Colorado field pilots. In: SPE Western Regional Meeting, Society of Petroleum Engineers, San Jose, California, 24–26 March.

Kibodeaux KR (2014) Evolution of porosity, permeability, and fluid saturations during thermal conversion of oil shale. In: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, Amsterdam, The Netherlands, 27–29 October.

Lee K, Moridis GJ and Ehlig-Economides CA (2016) A comprehensive simulation model of kerogen pyrolysis for the in-situ upgrading of oil shales. SPE Journal 21(5): 1612–1630.

Lee K (2014) Rigorous simulation model of kerogen pyrolysis for the in-situ upgrading of oil shales. Doctoral dissertation, Texas A & M University, College station, Texas, pp. 38–170.

10.

Moridis

(2014) TOUGH+ HYDRATE v1. 2 User’s manual: A code for the simulation of system behavior in hydrate-bearing geologic media, Berkeley: Lawrence Berkeley National Laboratory, pp. 9–74.

11.

Nimblett

Ruppel

(2003) Permeability evolution during the formation of gas hydrates in marine sediments. Journal of Geophysical Research: Solid Earth 108(B9): 1–17.

12.

Pruess

(1985) A practical method for modeling fluid and heat flow in fractured porous media. Society of Petroleum Engineers Journal 25(01): 14–26.

13.

Shen C (2009) Reservoir simulation study of an in-situ conversion pilot of Green-River oil shale. In: SPE Rocky Mountain Petroleum Technology Conference, Society of Petroleum Engineers, Denver, Colorado, 14–16 April.

14.

Symington W (2006) ExxonMobil’s electrofrac process for in situ oil shale conversion. In: 26th Oil Shale Symposium, Colorado School of Mines, Golden, Colorado, 16–20 October.

15.

Symington

Olgaard

Otten

(2010) ExxonMobil’s Electrofrac™ Process for In Situ Oil Shale Conversion, ACS Symposium Series, Oxford: Oxford University Press, pp. 185–216.

16.

Thoram S and Ehlig-Economides CA (2011) Heat transfer applications for the stimulated reservoir volume (SRV). In: SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, Denver, Colorado, 30 October–2 November.

17.

Vinegar H (2006) Shell’s In-situ conversion process. In: 26th Oil Shale Symposium, Colorado School of Mines, Golden, Colorado, 16–20 October.

18.

Wellington SL, Berchenko IE, De Rouffignac EP, et al. (2005) In situ thermal processing of an oil shale formation to produce a condensate. Google Patents (US 6880633 B2).