The shale samples analyzed in this study were collected from the Sichuan Basin, China. As illustrated in Figure 1, the Sichuan Basin, situated in the eastern Qinghai-Tibet Plateau and the northwestern Yangtze Plate, is a multi-cycle, rhombohedral sedimentary superimposed basin that has undergone multi-stage and multidirectional deep fault activity throughout its geological history (Liu et al., 2018). Over geological time, a series of marine shale formations were deposited in the Ordovician Wufeng Formation and Silurian Longmaxi Formation within the Sichuan Basin. The burial depth of the shale reservoir generally ranges from 2,000 to 5,500 m. The Wufeng Formation predominantly contains black graptolite-rich siliceous shales and ash-bearing siliceous shales throughout the region, while gray shales or silty shales are found near the paleo-uplift. The thickness of this shale sequence is generally less than 10 m (Nie et al., 2017). Conversely, the Longmaxi Formation is characterized by extensive deposits of black graptolite-rich, organic-rich shales with thicknesses ranging from approximately 85–105 m.

Table 1 summarizes the physical properties of the shale samples. Following established experimental protocols, porosity was measured using the helium method, and total organic carbon (TOC) was determined through high-temperature combustion of powdered samples. Specific surface area and pore size were determined through liquid nitrogen (N₂) adsorption experiments. Shale pore structure was characterized by analyzing the adsorption and desorption curves of N₂ within shale pores under low temperatures (77.4 K) and pressure. These data enable the calculation of specific surface area and pore size. Specific surface area was calculated using the BET equation (Brunauer et al., 1938), while pore size was determined through the BJH equation (Barrett et al., 1951). Whole-rock mineral composition and the relative content of clay minerals were analyzed using an X-ray diffractometer. Images of shale samples are provided in Figure 2.

Table 2 details the isothermal adsorption data and conditions gathered in this research. It includes data for samples No. 1-5 from Table 1, collected under varying temperatures and water saturation levels. This information was analyzed to develop models dependent on both temperature and moisture. Specifically, data from samples 1 and 2 were used to explore how temperature affects the adsorption mechanism of shale, while samples 1 and 3 provided insights into the impact of water saturation. Additionally, samples 4 and 5 were examined to understand the combined effects of temperature and moisture. Through comprehensive analysis of these samples, a model accounting for both temperature and moisture influences was developed.

The isothermal adsorption method is widely recognized as the standard approach for evaluating shale’s adsorption properties. This technique evaluates the gas adsorption capacity of shale per unit mass by maintaining constant temperature and varying pressure. In this study, the volumetric method was used to quantify shale’s adsorption capacity per unit mass. The volumetric technique involves placing the shale sample in a chamber, introducing helium to fill the non-adsorbed space, followed by methane gas. Methane is partially adsorbed by the shale, reaching equilibrium over time. The adsorption capacity per unit mass is determined by the change in gas volume within the chamber. At a constant temperature, methane is adsorbed at different pressures, forming a characteristic adsorption layer on the shale surface. However, as the adsorbed phase volume cannot be directly measured, the absolute adsorption amount remains inaccessible. The excess adsorption amount (or experimental adsorption amount) is measured, with the relationship between excess and absolute adsorption amounts provided by the following conversion formula: m a b s = m e x / 1 − ρ g / ρ a = m e x + ρ g V a

Where m _abs represents the absolute adsorption amount, g; m _ex is the excess adsorption amount, g; ρ _g is the density of methane, mol·m⁻³; ρ _ɑ refers to the density of the adsorbed phase, mol·m⁻³.

Rangarajan et al. (1995) applied the mean-field approximation to general density functional theory, leading to the formulation of the simplified local density (SLD) model. This model posits that interactions between adsorbate molecules, as well as those between adsorbates and adsorbent surfaces, collectively govern the adsorption process. Specifically, fluid-fluid interactions among adsorbate molecules and fluid-solid interactions between adsorbates and the adsorbent surface are both critical to the adsorption process. Interactions among adsorbate molecules are modeled using an equation of state (EOS), whereas those between adsorbates and pore walls are represented by potential energy functions. As shown in Figure 3, the SLD model simplifies the pores of the adsorbent into slit pores. The adsorbates with z from one pore wall are located between the pore walls with a pore width of L, and the adsorbates are simultaneously subjected to the common forces of both pore walls and other adsorbates.

At adsorption equilibrium, the chemical potential of adsorbates at a specific position z includes both fluid-fluid and fluid-solid chemical potentials, equating to the chemical potential in the bulk phase. μ z = μ f f z + μ f s z = μ b u l k (1)

Where μ(z) is the chemical potential at the z position in the pore, J·mol⁻¹; z is the distance between the adsorbates and the pore wall, nm; μ _ff (z) is the fluid-fluid interaction chemical potential at position z in the pore, J·mol⁻¹; μ _fs (z) is the fluid-solid interaction chemical potential at position z in the pore, J·mol⁻¹; μ _bulk is the bulk chemical potential in the pore, J·mol⁻¹.

According to thermodynamic equilibrium, chemical potential can be expressed by fugacity in nanopores: μ b u l k = μ 0 T + R T ⁡ ln f b u l k f 0 (2) μ f f z = μ 0 T + R T ⁡ ln f f f z f 0 (3)

Where μ ₀ (T) represents any reference state chemical potential, J·mol^-1; ƒ₀ refers to the fugacity of any reference state, Pa; ƒ_bulk and ƒ _ff (z) are respectively the bulk fugacity and the adsorption phase fugacity at z position in the pore, Pa.

In nanopores, the adsorbates are subjected to the force of the pore walls on both sides, and the chemical potential generated can be expressed as follows (Rangarajan et al., 1995): μ f s z = N A ψ f s z + ψ f s L − z (4)

Where N _A is Avogadro’s number; Ψ ^fs (z) and Ψ ^fs (L-z) are the potential energy generated by the interaction between adsorbates at position z in the pore and the pore walls on both sides, J.

Equations 1–4 can be obtained simultaneously: f f f z = f b u l k ⁡ exp − ψ f s z + ψ f s L − z k T (5)

Where k refers to Boltzmann’s constant, 1.3806505 × 10⁻²³ J·K⁻¹.

The potential energy generated by the interaction between adsorbates at position z in the pore and the pore walls on both sides can be calculated by Lee’s 10–4 potential model (Lee and Brenner, 1988): ψ f s z = 4 π ρ a t o m s ε f s σ f s 2 σ f s 10 5 z + σ s s / 2 10 − 1 2 ∑ i = 1 4 σ f s 4 z + σ s s / 2 + i − 1 σ s s 4 (6) where, ε f s = ε f f × ε s s (7) σ f s = σ f f + σ s s 2 (8)

Where ρ _atmos is the density of carbon atoms, 38.2 atoms·nm⁻²; ε _fs is the fluid-solid interaction energy, J; ε _fs is the fluid-fluid interaction energy, J; ε _ss is the solid-solid interaction energy, J; σ _fs is the diameter of fluid-solid, m; σ _ss is the interplanar distance of carbon, 3.35 × 10⁻¹⁰ m.

In the SLD model, the fluid equation of state is used to determine parameters like fluid fugacity, bulk fugacity, and fluid density. In this study, after evaluating various fluid equations of state, the PR-EOS was selected for computing these parameters. The equation is given by Subramanian et al. (1995): P R T ρ b = 1 1 − b ρ b − a T ρ b R T 1 + 2 b ρ b − b 2 ρ b 2 (9) where a T = 0.457535 α T R 2 T c 2 P c (10) b = 0.077796 R T c P c (11)

The term α(T) is expressed as (Gasem et al., 2001): α T = exp 2.0 + 0.8145 T r 1 − T r 0.134 + 0.508 ω − 0.0467 ω 2 (12) Where P represents pressure, Pa; ρ _b is the bulk density, mol·m⁻³; ɑ and b are respectively the attraction parameter and covolume, m³·mol⁻¹; ω is the eccentricity factor, T _c is the critical temperature, K; P _c is the critical pressure, Pa. The related physical properties of methane are shown in Table 3 Reid et al. (1987).

In the PR-EOS, the bulk fugacity is calculated as follows. ln f b u l k P = b ρ b 1 − b ρ b − a T ρ b R T 1 + 2 b ρ b − b 2 ρ b 2 − ln P R T ρ b − P b R T + a T 2 2 b R T ln 1 + 1 + 2 b ρ b 1 + 1 − 2 b ρ b (13)

Similarly, considering the fluid-fluid interaction, the fluid fugacity in situ of the nanopore is calculated as follows: ln f f f P = b ρ z 1 − b ρ z − a z ρ z R T 1 + 2 b ρ z − b 2 ρ 2 z − ln P R T ρ z − P b R T + a T 2 2 b R T ln 1 + 1 + 2 b ρ z 1 + 1 − 2 b ρ z (14) Where ρ(z) refers to the fluid density in situ of the nanopore, mol·m⁻³; ɑ(z) is the attraction parameter, J·m³·mol⁻². ɑ(z) is the position function of adsorbate molecules in nanopores, which can be calculated by the formula proposed by (Chen et al., 1997).

By combining Equations 5–14, density distribution of shale gas in nanopores is thus calculated. The expression of excess adsorption amount is as follows: n e x = 11.2 Z A ∫ l o w e r u p p e r ρ z − ρ b d z (15)

Where n _ex is the excess adsorption amount, m³·g⁻¹; A is BET specific surface area, m²·kg⁻¹; The lower limit and upper limit of the integral are respectively 3/8σ _ff and L-3/8σ _ff (Qi et al., 2019); Z is the gas compression factor.

To determine the in situ fluid density using the SLD model, the pores were divided into intervals of 0.01 nm. The local in situ fluid density for each interval was calculated, and the excess adsorption was determined using the composite Simpson’s method for integration.

The aim of this research is to develop a shale gas adsorption model that considers variations in temperature and moisture content, using data from isothermal adsorption experiments alongside the SLD model. This approach enables the direct calculation of shale gas adsorption under different conditions. The problem of parameter regression is reframed as a least-squares error minimization problem with two input and two output parameters (Golsanami et al., 2019; Yang et al., 2016). To address this, the paper uses the Levenberg-Marquardt (LM) algorithm to solve the least squares problem.

The core idea of the LM algorithm is its integration with the trust region method, allowing the objective function to satisfy the linearization condition within a specific interval, thereby effectively solving nonlinear least squares problems. The core idea is to balance the convergence speed and stability of the algorithm by dynamically adjusting the damping factor, combining the advantages of Gauss-Newton method and gradient descent method. The calculation formula of LM method is: S x = ∑ i = 1 m y i − f i x 2

Where y _i is the experimental result and f _i (x) is the predicted value of the model.

Temperature is a condition that cannot be ignored when studying the adsorption characteristics of shale. In this section, methane adsorption characteristics of two shales at different temperatures (303.15K-353.15K) are studied. The isothermal adsorption experiment results are shown in Figure 4 and Table 4. In this paper, the mean absolute error generated when calculating model regression by %AAD is used to characterize model accuracy: % A A D = 100 × A B S ∑ i = 1 N n i , c a l c u l a t e d − n i , t e s t / n i , t e s t N

Where n _i,calculated and n _i,test are respectively calculated adsorption capacity and experimental adsorption capacity at adsorption equilibrium, m³·g⁻¹. N represents the number of equilibrium adsorption points.

As can be seen from Figure 4, when the temperature remains unchanged, the excess adsorption amount rises sharply with the increase of pressure at the low-pressure stage. With the pressure increasing gradually, the excess adsorption amount increases slowly and reaches the maximum value at about 12 MPa. As the pressure continues to increase, the excess adsorption amount gradually decreases. The excess adsorption amount gradually decreases with the increase of temperature at the low-pressure stage. In the high-pressure phase, there is a reversal phenomenon, that is, with the increase of temperature, the excess adsorption amount gradually increases. This is because the adsorption of methane by shale is physical adsorption, which is an exothermic process. As the temperature rises, the kinetic energy of adsorbed methane increases, resulting in the reduction of the force between methane and the solid surface, thus weakening the adsorption capacity of shale (Chen et al., 2019; Heller and Zoback, 2014; Huang et al., 2022; Li et al., 2016). In addition, it is found that the adsorption of gas in shale is positively correlated with the number of micropores. The increase of temperature will change the pore structure of shale, resulting in the blocking or reduction of the number of micropores, and further weakening the methane adsorption capacity of shale (Li et al., 2021a; Mu et al., 2022).

Water generally exists in shale reservoirs, which can affect the gas adsorption capacity of shale reservoirs. The initial water saturation range of shale reservoir is 0%–70% (Ali et al., 2022; Ren et al., 2019). Therefore, this section studies the adsorption characteristics of shale in the 0%–70% water saturation range of two Longmaxi shales. The isothermal adsorption experiment data are shown in Figure 5. The SLD model is used to fit the data, and the regression results are shown in Table 4.

Figure 5 illustrates that with constant water saturation, the excess adsorption increases sharply as pressure rises during the low-pressure phase. As pressure continues to rise, the rate of excess adsorption increase progressively slows. After reaching its peak, the excess adsorption gradually decreases as pressure continues to increase. With constant pressure, the excess adsorption gradually decreases as water saturation increases. Unlike the effect of temperature on isothermal adsorption, inversion of excess adsorption does not occur during the high-pressure phase.

Water molecules inhibit methane adsorption on shale pore surfaces. The adsorption process of shale with water is shown in Figure 6. The interaction between liquid, gas and solid phases in shale pores affects the adsorption characteristics of methane molecules in shale (Li et al., 2016). At low water saturation, water and methane molecules compete for adsorption. However, the molecular force between water molecules and the pore surface of shale is hydrogen bond (Ma and Yu, 2022; Yang et al., 2020a), and the molecular force between methane molecules and the pore surface is van der Waals force. The hydrogen bonds are stronger than van der Waals forces. Water preferentially adsorbs and occupies the methane adsorption potential, which leads to the decrease of methane adsorption capacity. At high water saturation, water molecules form water films on the pore surface of shale and even produce capillary condensation in the nanopores, blocking the pore size. The production and thickening of water film narrow the diffusion path of methane molecules (Zhang and Yu, 2022), resulting in methane molecules unable to contact the pore surface of shale.

The effect of moisture on pore size of shale is another factor affecting the adsorption capacity of shale. Methane preferentially adsorbs in pores with strong adsorption capacity, and micro-pore filling occurs with the increase of pressure (Mu et al., 2022; Yang et al., 2020b) found that the maximum adsorption phase density occurred in the mesopore with a pore size of 4 nm. The competitive adsorption between water and methane is mainly concentrated in the mesopore (Zhang and Yu, 2022), and the affinity of pores to water is higher than that to methane in pores with pore size less than 30 nm (Ruppert et al., 2013; Sun, 2020) found this by studying pore widths in dry and moisture samples. Compared with dry samples, the content of pores with pore width less than 10 nm in moisture samples is significantly reduced, while the content of pores with pore width more than 10 nm is not significantly different, indicating that water in shale mainly affects pores with pore width less than 10 nm. Therefore, water content in shale decreases the mesopore content of methane, which leads to the decrease of adsorption capacity of shale.

Under actual reservoir conditions, high temperature and moisture coexist, which have synergistic negative effects on gas adsorption capacity. This section studies the methane adsorption characteristics of two shales under different temperature and water saturation states. The isothermal adsorption experiment results are shown in Figure 7. The SLD model and LM algorithm are used to fit the data (sample 5 lacks specific surface area and other experimental results, and is not fitted by the SLD model). The regression results are shown in Table 4. The relatively low %AADs listed in Table 4 indicate that the experimental results are in good agreement with those calculated using the SLD model.

According to Figure 7, during the change of methane pressure, the variation law of excess adsorption amount in shale is basically consistent with that described above. Compared with dry shale, the adsorption capacity of methane in moisture shale is significantly lower.

Temperature and moisture individually reduce the gas adsorption capacity of shale reservoirs; however, their coexistence results in combined effects. Previous studies have shown that water molecules preferentially occupy micropores in water-containing shale, while methane is adsorbed in micropores larger than 1.5 nm. Increasing temperatures can alter the pore structure of shale, leading to clogged or reduced pore sizes. In such cases, temperature and moisture may counteract each other as they primarily influence methane adsorption within the same shale pores. Gasparik et al. (2013) and Zou et al. (2019) also concluded that water and high temperature counteract each other based on their study of shale’s thermodynamic parameters.

According to the above analysis and discussion, the SLD model can well fit the adsorption capacity of shale under different specific conditions. According to Table 4, εss varies with water saturation and temperature, which indicates that there is a specific functional relationship between εss and them (Fitzgerald, 2005).found a linear relationship between ε_ss and temperature through extensive experiments: ε s s , T = ε s s , 0 + α T − T 0

Where ε _ss,T and ε _ss,0 are the solid-solid interaction energy at temperature T and T ₀ , respectively; α is the regression coefficient.

Studies have found that the relationship between excess adsorption amount and water saturation is not simply linear, but exponential decay (Chen et al., 2012; Li et al., 2021b).Therefore, to characterize the influence of moisture on ε_ss, the following exponential function is introduced: ε s s , w = ε s s , 0 ⁡ exp − δ S w

Where ε _ss,w and ε _ss,0 are the solid-solid interaction energy at temperature S _w and 0, respectively; δ is the regression coefficient.

Under actual reservoir conditions, high temperature and moisture coexist, which have synergistic negative effects on gas adsorption capacity. Therefore, the prediction model of shale gas content considering the influence of temperature, pressure and moisture is: ε s s , i = ε s s , 0 + α T − T 0 * ⁡ exp − δ S w (16)

Equation 15 and Equation 16 are used to globally fit all the excess adsorption data of different temperatures and water saturation in this study. The results of rock physics modeling using the LM algorithm are shown in Table 5 and Figure 8.

For the convenience of comparison, only the isothermal adsorption curves of shale at different temperatures and water saturations with corresponding experimental data are predicted. It can be seen from Table 5 and Figure 8 that the adsorption data calculated based on the temperature- and moisture-dependent model are very consistent with the experimental data, with %AAD between 1.8% and 8%, indicating that the model has a good fitting effect. The coefficients in the temperature- and moisture-dependent model shown in Table 5 may change due to the different geological conditions and reservoir properties of different shale gas reservoirs. This study aims to provide an idea for the establishment of a multi-factor shale gas adsorption model. According to this idea, the relationship between the adsorption capacity of shale in shale gas reservoir and various influencing factors is explored through experimental measurement, and a multi-factor adsorption model suitable for this reservoir is established, which can not only reduce the workload of experimental measurement on the premise of ensuring the calculation accuracy, but also provide guidance for the evaluation of adsorption gas.