Owing to their substantial reserves and the gradual exhaustion of conventional petroleum resources, in the past few decades there has been an interest in unconventional shale reservoirs as they represent large resources, cleaner energy option, and are seen as potential targets for CO₂ sequestration and thereby a potential bridge toward the transition to renewable energy resources (Jarvie et al. 2007; Godec et al. 2013; Heller and Zoback 2014; Hazra et al. 2018b, 2022; Du and Nojabaei 2019; Lyu et al. 2021; Chandra et al. 2022, 2023). Shale reservoirs commonly exhibit extremely low permeability (10^− 1 to 10^− 4 mD) (Wang et al. 2017; Vishal et al. 2019; Xu et al. 2023), making the extraction of shale gas from these reservoirs challenging through conventional methods. However, development of extraction technologies, specifically hydraulic fracturing (HF) has resulted in unprecedented development of shale gas exploitation in the past decade (Xu et al. 2023). Despite its substantial contribution to enhancing shale gas production, the hydraulic fracturing (HF) technique is not without drawbacks, including the potential for groundwater pollution and the risk of inducing earthquakes (Bao and Eaton 2016). Moreover, HF leads to extensive water penetration and pore space trapping (Liu et al. 2020). It is crucial to continuously explore novel methods for the extraction of shale gas. In-situ heating technique such as oil shale retorting is used as alternate to HF for extraction of shale oil/gas (Liu et al. 2020). The technique of in-situ oil shale retorting, dating back to the 1970s (Krumm 2014), involves drilling wells directly into the deposit from the ground. Additional heat, such as electrical heating, is applied within the formation to facilitate in-situ rectification of oil (Kang et al. 2020).

Several studies have explored methods to improve hydrocarbon recovery from gas shales through thermal treatments, including techniques such as in-situ combustion and pyrolysis (Chaprio and Bruining, 2015; Chen et al. 2016, 2018; Chandra et al. 2021, 2023; Hazra et al. 2023). Combustion or pyrolysis induces structural changes in kerogen, leading to the regeneration and reopening of pre-existing micro- and mesopores (Kang et al. 2020; Chandra et al. 2021, 2023; Hazra et al. 2023). This process enhances effective porosity, consequently improving hydrocarbon flow within the matrix. While the effects of pyrolysis on oil shales and its consequences for in-situ thermal treatment-based oil recovery have been widely studied (Tiwari et al. 2013; Sun et al. 2015; Bai et al. 2017; Lei et al. 2021; Taheri-Shakib and Kantzas 2021; Guo et al. 2022; Shi et al. 2023), there has been limited attention to its impact on gas shales and the consequent implications for in-situ thermal treatment-based hydrocarbon recovery.

Recently, a new subsurface thermal treatment technique, known as underground coal thermal treatment (UCTT) has emerged for stimulating deeply buried coals to produce lighter hydrocarbons (Krumm 2014; Kelly et al. 2016; Zhang et al. 2017; Shi et al. 2022). In this method, energy is applied to heat the coal seam through pyrolysis, and the resulting decomposition products are collected as gaseous and liquid fuels. A significant portion of the carbon remains below the surface, creating a carbon matrix, which can additionally be targeted for CO₂ sequestration (Krumm 2014; Shi et al. 2022). Compared to underground coal gasification (UCG), UCTT has a number of advantages due to the non-oxidizing environment. Additionally, pyrolytic semi-coke or the coal char produced by UCTT serves as good storage vehicle for CO₂ geological storage (Shi et al. 2022). Similar application on organic-rich shale horizons could also become a viable method for recovery of lighter hydrocarbons and creating an in-situ shale-char with higher surface areas, which can be conducive to CO₂ storage. Several laboratory scale experimental studies have shown the impact of heat treatment on shale petrophysical properties with implications for gas recovery and storage (Chen et al. 2017, 2018; Li et al. 2022; Chandra et al. 2021, 2023; Zhang et al. 2023; Zhuoke et al. 2023). Chen et al. (2017) conducted oxic heating of thermally matured shale samples, and observed removal of organic pores at temperatures of 400–500 °C. Chen et al. (2018) observed a reduction in surface area and an increase in the mean diameter of pores in shale samples subjected to oxic heating at 500 °C. Li et al. (2022) observed that with increasing thermal treatment temperatures, there was a continuous increase in the total pore volume, quantity, and ratio of mesopores and macropores in Longmaxi shale samples. Chandra et al. (2021) conducted both oxic and anoxic heating on Permian shale samples of Raniganj basin and witnessed a notable rise in surface area and pore volume during oxic heating, while anoxic heating led to a significant reduction in these parameters. Similar observations were reported by Chandra et al. (2023) while also highlighting the augmentation of surface area and pore volume through thermal treatment enhances the feasibility of long-term CO₂ storage in shale formations. Zhang et al. (2023) found that the Brunauer–Emmett–Teller specific surface area (BET SSA) of Bohai Bay basin shale samples initially decreased, then increased, and finally decreased again as temperature rose from 110 to 500 °C. In contrast, the average pore diameter demonstrated the opposite trend. Zhuoke et al. (2023) studied the effect of total organic carbon (TOC) on evolution of pore system in immature shale samples after heat treatment and observed that higher organic content correlates with increased hydrocarbon expulsion strength and pore volume growth. The pores exhibited pronounced fractal characteristics, with fractal dimension ranging between 2.397 and 2.636, indicating extremely small and complex pore structures with strong heterogeneity. Hazra et al. (2023) studied the impact of oxic heating (at 325, 425, and 525 °C) on the evolution of pore structures in the shale samples of contrasting thermal maturities. They observed that thermal treatment at higher temperatures significantly altered mesopore and micropore structures in shales.

In this study, two shale samples, with contrasting organic richness (TOC: 16.46 wt% and 28.98 wt%) collected from the depth intervals of 426-426.40 and 580–580.90 m, of a borehole drilled at Surni Coal Block of Rajmahal Coalfield, Jharkhand, India) belonging to the Lower Permian (Epoch- Cisularian; Stage: Kungurian) Barakar Formation from Rajmahal basin, India, are thermally treated (pyrolysis) at different temperatures to determine the temperature protocols needed for liberating hydrocarbons from the shales, and additionally explore the evolution of pore structural facets due to thermal treatment for CO₂ sequestration. In general, shales from Rajmahal are thermally immature (equivalent to sub-bituminous coal rank; Sethi et al. 2023) and hence presents the opportunity to thermally stimulate them to extract cleaner hydrocarbons and create sites for CO₂ sequestration. Examining thermally immature shale and its thermally treated counterparts also allow discerning the changes in the geochemical and pore properties as shale progresses towards maturity. Our results indicate that when the samples were heated at elevated temperatures for 30 min, the hydrocarbons are liberated from the kerogen, resulting in substantial increase in micropore surface area and volume of the samples, implying an increased CO₂ sequestration potential in thermally stimulated shales.

2.1 Sample types and thermal treatment

Two organic-rich shale samples belonging to the Lower Permian Barakar Formation from Rajmahal basin, India were selected for this study. The shales were controlled-crushed using pastel and mortar to achieve 212 microns size, which were then thermally treated, after removing a small aliquot representing the original sample. For thermal treatment, the samples were loaded in stainless steel containers, with openings in the lid, and treated in a tubular furnace, with nitrogen gas flow to maintain inert conditions. The samples were heated at temperatures of 350, 500 and 650 °C for 30 min. The raw samples (untreated aliquots) are denoted as X-RT and Y-RT, while the splits thermally treated at 350, 500 and 650 °C for 30 min are labelled as X-350/Y-350, X-500/Y-500, and X-650/Y-650, respectively.

2.2 Rock-eval

Rock-Eval 6, manufactured by Vinci Technologies was used to analyse the hydrocarbon generation properties of the raw shales and the thermally treated splits. A modified version of Rock-Eval’s “basic/bulk-rock method” was used to analyse the samples, with 750 °C set as the final oxidation temperature. This approach followed the protocols developed by Hazra et al. (2017, 2019) to ensure reliable data generation for Type III-IV kerogen rich samples.

2.3 Low pressure gas adsorption

Low pressure gas adsorption (LPGA) analysis of the raw shale samples (X-RT and Y-RT), and its thermally treated counterparts (X-350/Y-350, X-500/Y-500, and X-650/Y-650) was conducted for evaluation for their pore structure characteristics, including surface area, pore volume, and size distribution. Samples crushed to 212 microns were used for conducting the experiments. Additionally, prior to conducting the experiments, the samples were subjected to degassing for 3 h at 110 °C (Hazra et al. 2018a; Singh et al. 2021). Quantachrome® Autosorb iQ apparatus was used for the analysis with Nitrogen (N₂) and carbon dioxide (CO₂) used as adsorbates. The relative pressure of the adsorbate (P/P₀) was increased steadily over the course of LPGA investigation until it reached the condensation pressure that is unique to the adsorbate’s temperature. In this case, P₀ represents the adsorbate’s condensation pressure, and P represents the saturation pressure at each pressure point. The LPGA experiment using N₂ was carried out at 77 K temperature and 1 bar of saturation pressure. A comprehensive set of 40 data points was recorded across the pressure ratio (P/P₀) spectrum, ranging from 0.01 to 0.99. The experiments involving LPGA using CO₂ were executed under precisely controlled conditions within a water bath sustained at a consistent temperature of 273 K. The experiments included P/P₀ values between 0.0005 and 0.03. Nitrogen (N₂) isotherms were also utilised to assess the Barret-Joyner-Halenda (BJH), pore size distribution (PSD), Brunauer, Emmett and Teller (BET) specific surface area (SSA) and Frenkel-Halsey-Hill (FHH) fractal dimensions. The specific surface area was determined using the multipoint BET equation. Through linearization of the BET equation, the slope and intercept reveal monolayer adsorption behavior, facilitating precise calculations of surface area based on the nitrogen monolayer adsorption volume and molecular weight. The Barret-Joyner-Halenda (BJH) model was used to calculate pore size distribution from nitrogen desorption isotherms (Barrett et al. 1951). Utilizing the Kelvin equation (Pirngruber 2016), it relates pore size to relative pressure, deriving differential pore volume and providing insights into the distribution of pore sizes within the shale. The micropores surface area was determined from CO₂ adsorption isotherms using the Dubinin-Radushkevich (DA) equation (Dabrowski 2001). Graphical analysis such as plotting ln(Q) against ε², where Q represents adsorbed gas, and ε signifies Polanyi potential, provides the value of \(\beta\) (adsorption potential). The micropore surface (S_micro) area was then calculated applying the equation:

The Dubinin-Astakhov (DA) method was used to calculate the micropore volume from the equation (Günay et al. 2007):

Where, \({V}_{0}\) is a material constant, \(\beta\) is the adsorption potential, \(E\) is the characteristic energy of adsorption, and \(P/{P}_{0}\) relative pressure. The integration over the natural logarithm of the pressure ratio, derived from adsorption isotherm data yields micropore volume. The density functional theory (DFT) method was used to assess the micropore size distribution (Thommes et al. 2015). The foundation of the DFT model lies in utilizing molecular statistical thermodynamics to determine the specific adsorption quantity within a defined pore range at a given experimental temperature and pressure. This is achieved by solving the equation for grand thermodynamic potential, which is intricately linked to the distribution of gas density within a specific pore space.

The low pressure N₂ gas adsorption data was used to calculate fractal dimensions of the shales and their thermally treated counterparts following the Frenkel-Halsey-Hill (FHH) adsorption isotherm equation (reference). The FHH technique requires the construction of the lnN versus ln(-lnX) graphic using the N₂ adsorption isotherm data, where X denotes P/P₀. The slope of a straight line fitted to that data trend is used to calculate the fractal dimension (D) applying Eq. (3)

Where, a is a constant and (D-3) is the slope of the best-fit line.

Ideally, two specific fractal dimensions are typically distinguished: (1) D1 (surface fractal dimension, calculated at the lower relative pressure range where the nature of adsorption is contemplated to be governed by van der Waals forces); and, (2) D2 (pore structural fractal dimension calculated at higher relative pressure range where the nature of adsorption is governed by capillary condensation (Khalili et al. 2000; Yao et al. 2008).

2.4 Field emission scanning electron microscopy

The two shale samples (X-RT and Y-RT) and their thermally treated counterparts (X-350/Y-350, X-500/Y-500, and X-650/Y-650) of 212 microns size after heating were gold coated for a duration of 2 min before being introduced into the Carl Zeiss Merlin VP Compact field emission scanning electron microscope (FE-SEM). Sample analysis was conducted using both Secondary Electron (SE2) and Angle Selective Backscattered (AsB) modes. The primary objective of this analysis was to discern and comprehend the structural alterations within the thermally treated shale samples in direct comparison to their untreated counterparts.

3.1 Changes in thermal maturity with increasing temperatures of treatment

Table 1 presents the results of the two shales (X-RT and Y-RT) and their thermally treated counterparts. For the raw shales, X-RT and Y-RT, T_max of the samples can be observed to be 430 and 429 °C, respectively, indicating them to be thermally immature (Peters and Cassa 1994). For these two samples, the S2 was found to be 17.89 and 24.27 mg HC/g rock, respectively. With increasing temperatures of treatment, the S2 of the samples systematically decreased, which is caused due to expulsion of hydrocarbons from the samples. Figures 1 and 2 shows the S2 curves of the two shales and their thermally treated splits. While for the raw shales and the splits treated at 350 °C (Figs. 1a and 2a), the S2 curves (and T_max) were observed to be reliable, for the splits treated at 500 and 650 °C (Figs. 1b and 2b), the nature of the S2 curves can be seen to be highly irregular and consequently the T_max calculated from these S2 curves are unreliable. This phenomenon primarily occurs as a result of the systematic expulsion of greater amounts of hydrocarbons from the samples with increasing temperatures of treatment.

S2 pyrograms of sample X-RT and its thermally treated counterparts: a X-RT and X-350 b X-500 and X-650

S2 pyrograms of sample Y-RT and its thermally treated counterparts: a Y-RT and Y-350 b Y-500 and Y-650

For heat-altered and low-TOC shales, characterized by irregular S2 curves and consequently unreliable T_max, Hazra et al. (2021) have recently introduced Rock-Eval S4T_peak as proxy for thermal maturity. The S4T_peak is calculated from the temperature peak observed in the S4CO₂ curve, generated during the Rock-Eval oxidation stage. For both the shales, with increasing temperatures of thermal treatment, the S4T_peak was observed to increase (Fig. 3). This indicates that with increasing temperatures of thermal treatment, as hydrocarbons are expelled from the samples, the remaining organic matter within the shales become more and more aromatized, and thereby requiring higher temperatures to react/crack and generate CO₂/CO during oxidation. The results also further substantiate the usage of Rock-Eval S4T_peak as a thermal maturity proxy, especially applicable for low-TOC and heat-altered shales where T_max from pyrolysis stage is less reliable. It can be effectively summarized from the Rock-Eval results that with increasing temperatures of treatment, hydrocarbons are systematically removed from the samples signified by lowering of S2 and PC, and the resultant shale-chars successively requiring higher temperatures (i.e., higher S4T_peak) during oxidation to decompose.

S4CO₂ oxidation graphics of sample a X-RT and b Y-RT along with their thermally treated counterparts, including corresponding S4T_peak values at different thermal treatment temperatures

One of the main advantages often cited for in situ thermal treatment of coal seams (under inert conditions) in comparison to underground coal gasification is the potential for extracting only lighter hydrocarbons, leaving a significant amount of carbon as char and creation of microporous structures in situ, which can subsequently be used for sequestrating CO₂ (Shi et al. 2022). Rock-Eval results from the thermally treated splits of the shales clearly substantiates that with increasing temperatures of treatment, the pyrolyzable carbon (PC) is significantly reduced due to liberation and removal of hydrocarbons. On the other hand, the residual carbon (RC) of the thermally treated splits are substantially high (albeit some reduction with increasing temperatures), indicating that significant portion of the carbon exists within the shale-char created from thermal treatment. The increasing S4T_peak from the oxidation stage of Rock-Eval also indicates the aromatization of the samples and formation of successively more aromatized shale-char owing to thermal treatment. The shale-char generated in situ through thermal treatment, may thus offer additional sites for storage of injected CO₂.

3.2 Inferences from low pressure N₂ gas adsorption

3.2.1 Pore structural parameters

Figures 4 and 5 plots the low-pressure nitrogen gas adsorption isotherms of the raw shale samples, and its thermally treated counterparts. For the two raw shales X-RT and Y-RT used in this study, the nitrogen gas adsorption derived BET SSA was observed to be 18.46 and 14.57 m²/g, with average pore radius being 8.26 and 7.78 nm, respectively.

Adsorption and desorption isotherms of a X-RT shale and its thermally treated counterparts: b X-350, c X-500, and d X-650, obtained from LPGA analysis using nitrogen gas

Adsorption and desorption isotherms of a Y-RT shale and its thermally treated counterparts: b Y-350, c Y-500, and d Y-650, obtained from LPGA analysis using nitrogen gas

The shape of the low pressure nitrogen gas adsorption isotherm reflects the type of pore present within the samples (Kuila and Prasad, 2013). The IUPAC classifies adsorption isotherms into 6 types and a detailed description of the classification of isotherms can be found in Sing (1985) and Rouquerol et al. (1998). For both the samples (Figs. 4 and 5), significant hysteresis was observed, implying occurrence of capillary condensation within the mesopores (Sing, 1985). However, none of the isotherms show any development of a plateau at the highest relative pressures, which is typical of type IV isotherms for mesoporous materials, with the plateau implying completion of mesopore filling (Kuila and Prasad, 2013; Hazra et al. 2018a). Rather the isotherms show development of steep slopes at the highest relative pressures, implying presence of macropores (Kuila and Prasad, 2013; Hazra et al. 2018a). Rouquerol et al. (1998) designated such isotherms as Type IIB i.e., samples marked by presence of both mesopores (causing hysteresis) and macropores (causing absence of plateau like those shown by Type IV materials). Consequently, the adsorption isotherms for the studied samples are identified as Type IIb (Figs. 4 and 5). The shape of the pores is inferred indirectly through the analysis of the N₂ adsorption isotherm, particularly focusing on the slope and curvature of the adsorption branches. The hysteresis loops in the desorption branch can indicate the pore shapes. De Boer (1958) classifies hysteresis loops into five types, each correlated with distinct pore shapes. Specifically, Type A hysteresis is linked to cylindrical pores, Type B to slit-shaped pores, Types C and D hysteresis are indicative of wedge-shaped pores, and Type E associated with bottle neck pores (Labani et al. 2013). Both shale samples exhibit Type B hysteresis loops (Figs. 4 and 5), suggesting the prevalence of predominantly slit-shaped pores.

As the thermal treatment temperatures increased, marked changes in the adsorption isotherm characteristics were observed for both the samples. Compared to the raw samples, the X-350 and Y-350 splits exhibited a significant reduction in adsorption capacity at the lowest experimental pressures. A similar reduction in adsorption capacity was also observed for the X-500 and Y-500 splits at the lowest experimental pressures. What this essentially indicates is that due to thermal treatment at 350 and 500 ºC, as the organic matter within the shales react and generate hydrocarbons, some structural reworking within the samples take place, which limits the entry of nitrogen in assessing the finer mesopore structures. On the other hand, for the sample splits thermally treated at 350 and 500 ºC, tightening of the hysteresis loops was noted, which indicates that the connectivity of coarser mesopores is improved during hydrocarbon generation due to thermal treatment. Consequently, a decrease in BET SSA and an increase in the average pore radius was observed for splits treated at 350 and 500 °C (Table 2).

The restricted entry of N2 into the pore spaces of the splits treated at 350 and 500 ºC is also manifested in the BJH pore size distribution (PSD) plot. Notably, for the X-350, X-500, Y-350 and Y-500 splits, the volume of pores in the smaller sizes (∼10 nm) is several times smaller compared to their respective raw splits (Fig. 6). Interestingly, compared to the raw sample and the other splits, for the X-650 split the BET SSA was observed to sharply increase (38.38 m²/g; even higher than X-RT), with the average pore radius being sharply reduced (3.25 nm; even lower than that of X-RT). However, for sample Y-650, BET SSA (4.84 m²/g) was observed to only marginally increase compared to Y-500 (4.61 m²/g), and thereby being much lower than Y-RT. Similarly, for the Y-650 split average pore radius was observed to be marginally lower compared to the Y-500 and Y-350 splits, but much larger than the Y-RT split. Careful consideration of the Rock-Eval data also shows contrasting behaviour of the two shales, despite of them having similar HI and thermal maturity levels (Table 1). S1 was first observed to increase from 0.10 mg HC/g rock (for X-RT) to 0.14 mg HC/g rock (for split X-350), followed by a sharp increase for the X-500 split (0.58 mg HC/g rock), followed by a sharp decrease for the X-650 split (0.14 mg HC/g rock). The sharp increase in S1 for the X-500 split could be due to generation of greater amounts of free hydrocarbons from the organic matter within the sample caused due to higher temperature thermal treatment. With further thermal stimulation (650 ºC), substantial portion of the free hydrocarbons are expelled from the sample. In contrast for sample Y-RT and the other splits, no sharp increase or decrease in S1 was noted, with near similar S1 for the Y-350, Y-500 and Y-650 splits but higher than the Y-RT split. While the BJH PSD plot shows opening of larger number of smaller mesopores in X-650, relative to X-350 and X-500 splits, the plot doesn’t fully substantiate the increase in BET SSA and decreased average pore size in X-650 relative to X-RT (Fig. 6a). However, DFT model showed larger volume of smaller pores in X-650, relative to X-RT. Figures 6c and d plot the comparative PSD of X-RT and X-650 splits, derived using DFT model. DFT method, unlike BJH adsorption model, takes into account the effect of adsorption potential on location of pore condensation transition and thereby is more suitable for characterizing narrow mesopores using nitrogen adsorption, which is generally underestimated in the BJH adsorption model (Chandra and Vishal 2020; Pang et al. 2021). In Fig. 6c, it can be clearly seen that the volume of pores [dV/dlog(D)] till 4 nm width is larger than in X-650, relative to X-RT. Similarly, the cumulative volume of pores at the narrow mesopores ranges in X-650 is manifolds higher than that is X-RT (Fig. 6d), although the overall pore volume is lower. This substantiates the higher BET SSA and lower pore size of X-650, relative to X-RT.

PSD plots for sample a X-RT and its thermally treated counterparts, and for b Y-RT and its thermally treated counterparts derived using BJH adsorption model c PSD plots of sample X-RT and d its thermally treated X-650 split derived using the DFT model

3.2.2 Changes in fractal dimensions due to thermal treatment

Frenkel − Halsey − Hill (FHH) fractal model was used for determining the fractal dimensions of the shales and to understand the evolution of pore complexities as thermal treatment temperatures increased. The details and how the FHH model is used to calculate fractal dimensions (D) are presented in Yang et al. (2014) and Li et al. (2016). In general, fractal dimensions are recognized to be of two types- D1 and D2. Smaller mesopores (2 to 8 nm diameter) gets saturated inside the relative pressure stretch of 0.01 to 0.5 and the behaviour here is controlled by the Van der Waals forces; fractal dimensions calculated in this domain is represented as D1. On the other hand, larger mesopores (> 8 nm to 50 nm diameter) are saturated in the relative pressure range of 0.5 to 0.99 through the mechanism of capillary condensation; fractal dimensions calculated in this domain is represented as D1. For both the set of shale samples used in this study, significant alterations in fractal dimensions was noted with changing thermal treatment temperatures (Table 2).

The FHH plots of the two shales and the thermally treated counterparts are presented in Figs. 7 and 8. For sample X-RT and its splits, in general, fractal dimension D2 was observed to show better correlation with other pore properties (viz. BET SSA and average pore diameter; Fig. 9A and B). D1 and D2 first exhibited a decline from 2.40 to 2.75, respectively (in the case of X-RT), to 2.24 and 2.67 (for X-350). While D2 was observed to fall further for X-500 (2.61), D1 showed substantial increase (2.44). For split X-650, similar to other pore properties, both D1 and D2 showed sharp increase (2.78 and 2.93, respectively), indicating formation of complex porous structures owing to thermal treatment. For the other set of samples (Y series), the changes in fractal dimensions with changing temperatures of treatment were observed to be non-systematic and slightly complex, with D2 showing good correlation with other pore parameters and D1 showing little or no correlation (Figs. 9c and d). Y-RT showed the highest fractal dimensions indicating presence of complex porous structures. A reduction in both D1 and D2 was noted for Y-350, followed by a subsequent decrease in D2 and an increase in D1 for Y-500. For the Y-650 split, while D2 was observed to increase, D1 was observed to be lower than Y-500. The non-systematic variations in the evolution of fractal dimensions of both sets of shale samples used in this study, indicates certain complexities in understanding their pore structural attributes. This could be attributed to closure of preexisting pores with increased temperatures of treatment or restricted entry of nitrogen into the smaller pores, thereby producing non-systematic variations. It is also important to consider that, low pressure nitrogen gas adsorption presents certain challenges, especially when considered for evaluating pore structures in organic-rich rocks, especially for those with TOC > 10 wt% (Hazra et al. 2020). Previous research works focused on coals by several (viz. Gan et al. (1972), Parkash and Chakrabartty (1986), Mahajan (1984, 1991), Mastalerz et al. (2009) and Okolo et al. (2015)) have revealed very low N₂ gas adsorption derived SSA of coals, while the same samples when analysed using CO₂ as the adsorbent revealed very high micropore surface areas. It may be argued that the data from N₂ and CO₂ should not be compared as essentially different pore structures (mesopores vs. micropores) are assessed by the two adsorbents, however, the decreased N₂ adsorption derived surface areas for coals and organic-rich rocks (TOC > 10 wt%) is somewhat counter-intuitive (Hazra et al. 2020). Moreover, N₂ due to its low kinetic energy at -196 °C experimental temperatures present challenges in accessing/penetrating intricate pore structures in coals and also due to potential shrinkage of pores at such low experimental temperatures (Anderson et al. 1962; Walker et al. 1988; Mahajan 1991). Additionally, coaly organic matter is also known to have molecular sieve properties with large number of pores being lower than equal to 4 Å, which limits entry of N₂ at -196 °C, but allows passage of CO₂ (Lamond 1962; Thomas and Damberger 1976; Mahajan and Walker 1978). The shales used in this study are marked by higher TOC content (16.46 wt% and 28.98 wt%), and are collected from the coal-bearing Barakar Formation of Rajmahal basin. Previous study in this area indicates the predominance of type III-IV kerogen within the shales (Sethi et al. 2023). For the studied thermally stimulated organic-rich shales, it is entirely possible that certain sieving effect is imparted by the coaly organic matter, resulting in reduction in surface areas with increasing temperatures of treatment. Similarly, the variable and non-systematic variations of fractal dimensions, especially D1 (representing the smaller mesopores (2 to 8 nm diameter) within the relative pressure stretch of 0.01 to 0.5) also indicates certain complexities, which might be related to the lack of assessment of porous structures formed due to thermal stimulation by the N₂ at the experimental conditions. Comparison of the data produced from N₂ with that of CO₂ would indicate whether the non-uniform variation generated are related to the pore properties of the shales or due to limitations of N₂ in accessing some of the pore structures.

Fractal plots of sample a X-RT shale and its thermally treated counterparts: b X-350, c X-500, and d X-650

Fractal plots of sample a Y-RT shale and its thermally treated counterparts: b Y-350, c Y-500, and d Y-650

Cross-plots illustrating relationship between fractal dimensions, D1 and D2 with other pore properties for the X-series (a and b) and Y-series (c and d) of samples. Note: The increasing sizes of the sample symbols corresponds to increasing thermal treatment temperatures

3.3 Inferences from low pressure CO₂ gas adsorption

Low-pressure CO₂ adsorption experiments were conducted on both the raw shales and the thermally stimulated counterparts to investigate the changes in micropore properties of the samples, with changing temperatures of treatment. The findings from the low-pressure CO₂ adsorption are presented in Table 3. The isotherms, overall, can be categorized as Type I, indicating a greater adsorption rate at lower relative pressures. Sample Y-RT, marked by higher TOC content, was observed to show higher adsorption capacity, surface area, pore volume and reduced pore size in comparison to X-RT. This observation suggests the existence of larger volumes of micropores within the organic-rich sample. This is opposite to the observation from N₂ adsorption, where higher surface area and pore volume was observed in sample X-RT, relative to Y-RT. When both the samples were thermally treated at 350 °C, some drop in CO₂ uptake capacity was noted (Fig. 10). Although a similar drop in N₂ adsorption derived mesopore properties was also noted, however, the decrease from CO₂ adsorption was observed to be less sharp. Nevertheless, this indicates closure of some of the micropores that were present within the raw samples due to thermal treatment at 350 °C. Micropore size distributions calculated using the CO₂-based DFT model further corroborated these findings (Fig. 11).Compared to X-RT for the X-350 split, till ~ 0.45 nm pore width, new additional openings can be noted; however, beyond that the volume of pores in the X-350 split can be observed to be substantially lower than the X-RT split (Fig. 11). For Y-350, the volume of pores was noted to be lower than that of Y-RT, thereby indicating closure of some pores (Fig. 11).

Low-pressure CO₂ gas adsorption isotherms of the a X series and b Y series samples used in this study

Micropore size distributions for the X series (a and b) and Y series (c and d) samples calculated using the CO₂-based DFT model

Opposite to what was observed from N₂ adsorption, from CO₂ adsorption a sharp increase in CO₂ uptake capacity was observed for X-500 (~ 1.5 times higher than X-RT) and Y-500 (~ 1.4 times higher than Y-RT). Similarly, for these splits, micropore surface areas was observed to sharply increase (higher than the raw splits). This increase on CO₂ micropore surface area, CO₂ uptake capacity for the X-500 and Y-500 splits, indicates significant formation of micropores owing to thermal treatment at 500 °C, which is also corroborated from the micropore size distributions (Fig. 11). For the splits treated at 650 °C, the micropore surface areas and micropore volumes were observed to further increase, indicating formation of newer micropores due to thermal treatment at 650 °C. Micropore size distributions plots also revealed formation of greater volume of pores in X-650 and R-650, compared to X-500 and Y-500. Overall, the results from CO₂ adsorption in a way contradicts the findings from N₂ gas adsorption which may be indicative of the limitations of N₂ in accessing complex pore structures formed in the organic-rich samples owing to thermal treatment.

3.4 Implications from electron microscopy

Generally, observations from Scanning Electron Microscopy or Imaging do not actually validate or substantiate the pore properties derived from low pressure gas adsorption. The target of conducting Scanning Electron Microscopy of the thermally treated shales was not to correlate the observations from SEM with low pressure gas adsorption, but to independently examine the impact of thermal treatment on SEM observable porosity. It is fairly well established that as thermal maturity of organic matter in sedimentary rocks increase, microporosity (< 2 nm) increases and the importance of macroporosity is diminished (Okolo et al. 2015; Mastalerz et al., 2018). Previous research has documented that SEM observations can only probe/image pores with diameters ≥ 2 nm (i.e., meso- and macro-pores) (Liu et al., 2017; Mastalerz et al., 2018), and micropores (< 2 nm) which plays the most important role is adsorption of gases (both for unconventional gas production and CO₂ sequestration) are not identifiable under SEM. Consequently, comparison/correlation of porosity determined from SEM with that of other methods such as low pressure gas adsorption is largely uncertain and likewise the legitimacy of any calculation of storage capacity based on SEM-identifiable porosity is largely doubtful (İnan et al., 2018).

The SEM images of the two shale samples and their thermally treated counterparts show interesting evolution in pore structure. X-RT shows microcracks in OM matrix (Figs. 12a, b). Pore size in X-RT have a modal average of 182 nm with a maximum pore size of 2 μm. Pores are elongated with a high aspect ratio and are primarily hosted within the organic matter. Both pores and microcracks serve as gas adsorption sites. X-350 shows elongated lenticular pores similar to X-RT (Figs. 12c, d), however the opening of the pores are wider and has lower aspect ratio compared to pores in X-RT, suggesting that thermal treatment opens up the pre-existing pores. After further thermal treatment at 500 °C, the X-500 shale shows more circular pores (Figs. 12e, f), i.e., pores with lower aspect ratio as compared to X-RT and X-350. No newly formed pores were found within the scanned region, however the increase in width of the existing pores results in increase in volume of the existing pores. After further heating at 650 °C, the X-650 shale shows interesting evidence of pore collapse (Figs. 12g, h). Pre-existing pores separated by thin layer of organic matter gets ruptured and two or more pores merge to form a bigger pore. The orange circle in Fig. 12h indicate evidences of rupture of pore wall between multiple pores, which was not seen in lesser thermally treated shales. Aside from merging of pores, some distortion of the organic matrix is observed in some cases, indicating possibility of new pore formation or dilation of extremely small pores; however, the resolution of SEM images was unable to conclusively indicate that.

SEM images of shale X (a, b)- X-RT untreated (c, d)- X-350 after 350 °C thermal treatment (e, f)- X-500 after 500 °C thermal treatment and (g, h)- X-650 after 650 °C thermal treatment. The marked area (h) represents signatures of pores merging due to thermal treatment

SEM images of shale Y (a, b)- X-RT untreated (c, d)- X-350 after 350 °C thermal treatment (e, f)- X-500 after 500 °C thermal treatment and (g, h)- X-650 after 650 °C thermal treatment

Similar to X-RT, Y-RT also shows numerous microcracks and elongated pores in the organic matrix and clay minerals (Figs. 13a, b). There are abundant mesopores with a modal average width of 300 nm based on the SEM images of different magnification taken at various spots of the visible surface. Similar to shale X, shale Y also shows elongated pores with higher aspect ratio at lower temperatures which starts dilating along the short axis with increasing temperature (Figs. 13c, d). Unlike X, the pore dilation between 350 °C (Y-350) (Figs. 13c, d) and 500 °C (Y-500) (Figs. 13e, f) is apparently insignificant. Y-650 (Figs. 13g, h) shows circular pores with smoother pore walls as seen in X-650, however there was no evidence of pore wall rupture as seen in X-650. For both X and Y, multiple images demonstrated evidence of clay mineral pores, however these pore systems didn’t show significant change with increasing temperature. A summary of the pore size, aspect ratio and their change with temperature for both the shale samples are shown in Table 4.

The results essentially indicate that owing to thermal treatment (under inert conditions), as more and more hydrocarbons are expelled from the samples, the remaining organic matter becomes more aromatized (as revealed from the Rock-Eval results), and is accompanied by formation of more micropores (as revealed from CO₂ adsorption experiments). These microporous sites, can be used for sequestrating atmospheric CO₂. Consequently, similar to underground coal thermal treatment (Zhang et al. 2017; Shi et al. 2022), if deep underground organic-rich shales are thermally stimulated, it can result is generation of cleaner hydrocarbons and also create sites for CO₂ sequestration. While, the temperature protocols used in the experiments for this study indicate that at 650 °C almost all pyrolyzable component from the kerogen is generated and the shale-char created at this temperature has maximum micropore capacity for CO₂ storage, however, depending upon the thermal maturity stage of a specific shale horizon and the type of kerogen it consists of, the temperature of in-situ pyrolysis can vary, as thermally immature shales with greater concentration of reactive kerogen start generating hydrocarbons at lower temperatures (Singh et al. 2021; Hazra et al. 2023). On the other hand, shales marked by higher thermal maturities and presence of less reactive kerogen would need higher temperatures for reacting, generating hydrocarbons and creating microporous sites for CO₂ sequestration. Thus, establishing the heating programs for in situ shale conversion should be based on the kerogen type present and the thermal maturity stage of the horizon.

In addition to alterations in shale pore structure resulting from thermal treatment, the project outcomes may also be impacted by the thermal treatment-induced propagation of fractures. The distribution and interconnection of fracture pathways within shale formations play a crucial role in the preservation, concentration, and productivity of pyrolysis byproducts (Geng et al. 2017; Li et al. 2022). When subjected to high-temperature treatment, the interconnected nanopores within the organic matter undergo enlargement, giving rise to numerous micropores and fractures, ultimately leading to increased complexity in the shale’s fracture networks and a significant enhancement in the quality of the shale matrix for gas recovery and storage (Kang et al. 2016; Lei et al. 2021). Moreover, the crack pattern and thermal cracking efficiency can be notably influenced by the thermal expansion coefficients of various minerals as the thermal treatment temperature rises (Liu et al. 2020). With increasing temperatures of thermal treatment, the number and maximum aperture of fractures in shales increase (Geng et al. 2017). The fracture toughness of shale diminishes with rising in-situ temperature, while simultaneously, fracture toughness anisotropy increases; notably, fracture paths aligned with bedding orientations of 0°, 30°, and 90° exhibit minimal sensitivity to temperature variations, whereas those with 45° and 60° bedding orientations are affected (Guo et al. 2023). In general, the heat induced fracture development in shale result in reduced strength and elastic properties, (Vishal et al. 2022; Zheng et al. 2023). Thus, a number of other factors viz. geomechanical stability of the thermally stimulated horizon, propagation of fractures, impact of supercritical phase transition of sequestered CO₂ (SC-CO₂) and its interaction with thermally stimulated shales etc., can influence the outcomes for such projects, and future studies directed in these areas can help in targeting underground organic-rich shales for extraction of cleaner hydrocarbons and sequestration of atmospheric CO₂.

1. (1)

   The results from low pressure N₂ gas adsorption experiments were observed to be in contradiction to the results when CO₂ was used as the adsorbent. This could be attributed to the relative inability of N₂ to access complex structures formed in the thermally treated splits.

2. (2)

   With rising thermal treatment temperatures, micropore surface area and volume, were first noted to decrease (by ∼ 7% and 13%, at 350 °C), followed by sharp increase (by ∼ 30% to 35% at 500 °C and ∼ 41% to 63% at 650 °C), for both the samples, compared to the raw samples, indicating formation of new micropores and thereby increasing the sequestration potential.

3. (3)

   Results from Rock-Eval indicate that for the shales used in this study, at 650 °C almost the entire pyrolyzable component of the kerogen is liberated. Additionally, the maximum micropore capacity determined from CO₂ gas adsorption experiments for these splits indicate 650 °C to be the best temperature for in situ conversion and CO₂ sequestration. However, other factors viz. geomechanical stability of the treated horizon, SC-CO₂-shale interaction should be considered for such projects.

4. (4)

   The FESEM imaging of the two shales and their heated counterparts reveal that the average pore width and pore aspect ratio consistently increases and decreases respectively and with temperature. When comparing untreated and 650 °C heated counterparts, sample X and Y show 11x and 29x increase in pore width respectively. This is accompanied by a reduction in aspect ratio by 8x and 10x for sample X and Y respectively. This indicates that the newly formed pores are not only bigger, but also appear more circular on visible surface.