A method for simultaneously determining axial permeability and transverse permeability of tight reservoir cores by canister degassing test

Permeability anisotropy of shale and tight gas reservoirs is critical for applications in unconventional gas recovery, but laboratory measurements are still limited. This paper presents an experimental method for determining permeability anisotropy of shale and tight reservoirs. The method uses gas production data from a canister and an analytical solution of continuity equation in anisotropic core sample. Axial permeability and transverse permeability of the core are estimated by matching experimental data with analytical solution. The proposed method is verified by comparing with true values and calculated values through numerical simulations that cover variations in rock permeability. The method is applied to real data measured in canister degassing tests (CDT) involving two different directional shale core samples. Then, the experimental results are compared with the results from pulse pressure decay (PPD) method. Both the verification from numerically calculated results and the comparison between PPD results and CDT results exhibit the practicality of the proposed method. At last, the effects of permeability anisotropy, porosity, initial gas pressure, rock dimensions, and adsorption on the profiles of gas production are investigated.


| 1221
ZHAO et Al. method, [3][4][5][6] oscillating pulse method, [7][8][9] canister degassing test method, 1,10,11 and the Gas Research Institute (GRI) method. 12,13 These methods have several advantages: (a) short experimental duration and (b) the pressure signal or the cumulative flow as opposed to instantaneous flow rate is measured. With these proposed methods, the effects of stress, 14 fracture properties, 15 water content, 16 and gas adsorption 17 on shale permeability are investigated. However, since these permeability measurement methods are developed by assuming one-dimension flow, only axial permeability of a core sample or spherical permeability of crushed samples can be obtained. Therefore, most of these studies are limited to homogenous rock. It is, in fact, important to develop a method to measure the anisotropic permeability because shale and tight gas reservoirs are indeed anisotropic rocks. The permeability anisotropy, defined as the ratio between permeability parallel to bedding and perpendicular to bedding, is generally greater than one, but varies greatly from 1.2 to 1864. 4 for Longmaxi and Wufeng shales from the Sichuan Basin in South China. 18 The large variation in permeability anisotropy can be attributed to the composition, grain size, porosity, effective stress, and microfractures. [19][20][21] Mokhtari and Tutuncu 22 tested directional (vertical, 45-degree, and horizontal) Mancos shale in Colorado using unsteady-state techniques and showed that the presence of lamination and fractures caused significant directional dependency of permeability.
Due to the lack of appropriate method to measure permeability anisotropy, the most used method is to drill two test cores parallel to bedding and perpendicular to bedding, respectively, to measure their permeabilities using unsteady-state techniques. However, the artificial microfractures of different core samples cause this method to be unreliable. A huge error may exist because the permeability measurement in different directions cannot be based on the same core sample. 23 In this paper, an experimental method for simultaneously determining axial permeability and transverse permeability of tight reservoir cores is developed. The experiments, known as the canister degassing test, are originally performed for calculation of lost gas by Hosseini et al 24 and are further improved by Alfi et al 11 for permeability measurement on homogenous rock. We modified their work to simultaneously determine axial permeability and transverse permeability of core samples. The method is verified through numerical simulations and is successfully applied to Longmaxi shale core to estimate the axial permeability and transverse permeability. At last, the influence factors and limitations of the proposed method are discussed.

| METHOD
The core samples are saturated with gas at initial gas pressure first and then degasified inside a canister at a lower constant pressure. Gas production data are recorded by the flowmeter for calculation of axial permeability and transverse permeability.
All sides of the core sample (top, bottom, and radial surrounding surfaces) are open to flow at constant pressure. Figure 1 shows the schematic model of canister degassing tests. The apparatus is composed of a canister, a high-resolution flowmeter (made from OMEGA company, USA, which can measure gas flow rate from 0 to 10 ml per min, with precision of 0.01 ml), pressure transducer (made from OMEGA company, which can measure gas pressure of 0 to 100 MPa, with precision of 0.001 MPa), back pressure regulator, gas supplement system (consists of gas cylinder, air compressor, and gas chamber, which can provide gas pressure of 0 to 20 MPa), and vacuum system.
The procedure of canister degassing tests is summarized as follows: (a) preparation, place the sample into the canister, and vacuumize the system for 24 hr. (b) Gas saturation, open valve 2 to bring the canister a required pressure. The core is saturated with Helium to achieve the desired equilibrium pressure. (c) Gas drawdown, a lower constant pressure (or atmospheric pressure), is applied at the outer boundary to let the core gas release through a back pressure regulator. More details of the canister degassing tests can be found in Ref. [1] and Ref. [11].
Two core samples, drilled from an outcrop of the lower Silurian Longmaxi formation in the Sichuan Basin, China, were used for permeability tests. Figure 2 shows the prepared samples, which are 50 mm in length and 50 mm in diameter. To investigate the permeability anisotropy, one sample (namely vertical sample) is drilled parallel to the bedding, and another (namely horizontal sample) is perpendicular to bedding. The mineral compositions of the samples are quartz, potassium feldspar, albite, illite, chlorite, Calcite, and Pyrite, The total organic carbon (TOC) of the samples is 2.88%, and vitrinite reflectance (Ro) is 3.72%.

| Analytical solution
The model is developed based on the single-phase gas flow in a radially and vertically finite core sample, considering the effect of transient flow with pseudopressures. 24 For the case of a heterogeneous cylinder core, the mass-continuity equation can be written as follows: The sample is saturated with gas at initial gas pressure first and then degasified inside a canister at a lower constant pressure. The cumulative volume of degassed gas is chronologically recorded by the flowmeter for the permeability measurement where ϕ is porosity, ρ g is gas density, t is time, and v r and v x refer to the flow velocities in r direction and x direction, respectively. Note that gas adsorption is not considered in this mass-continuity equation.
To account for changes in the compressibility and viscosity of gas caused by pressure decline, the pseudopressure is used as 25 Based on the equation of state, the gas density is as follows: where M is molecular weight of gas, R 0 universal gas constant (8.314 J K -1 mol -1 ), and T is experimental temperature.
Assuming that Darcy's law prevails, gives where k r and k x are the permeabilities in r direction and x direction, respectively.
Replacing these parameters give The boundary and initial conditions are as follows: where L and R denote length and radius of the core sample. m p0 is the pseudopressure at the equivalent pressure, and m pi is the initial gas pseudopressure in the core sample. For convenience, a set of dimensionless parameters are defined as Replacing for the dimensionless parameters yield: where ν 2 is defined as The solution of Equation (9) with corresponding boundary and initial conditions for dimensionless pressure is given by 24 (1) where λ n is defined as the roots of λ n are of the zero-order, first-kind Bessel function: λ 1 = 2.40483, λ 2 = 5.52008, λ 3 = 8.65373, λ 4 = 11.7915, and λ 5 = 14.9309 and F(λ n ) is.
With the pressure profiles, gas produced from the sides of the core and the top-bottom portion can be calculated as follows: where q s and q tb are the gas flow from the sides of the core and the top-bottom portion, respectively. Adding these terms gives the total gas production from the core sample as follows: where V c is the volume of the core, η is the inverse of diffusivity coefficient, ω m is the Fourier transform eigenvalues, λ n Hankeltransform eigenvalue, and χ is the ratio of axial permeability to transverse permeability.
The gas production is analyzed by history matching the analytical solution to the measured test data. In the history-matching process, there are only two unknowns: transverse permeability and permeability ratio. All other parameters, including core dimensions and gas properties, are assumed to be known. Generally, parameter estimation problems can be transformed into least-squares problems. Many optimization methods are specially designed to solve least-squares problems developed from Newton's method. One very popular local optimization method is the Levenberg-Marquardt (LM) method. However, the LM method alone is not sufficient to solve our problem because Equation (15) is nonconvex problem, which means that it has many local minima and that the global minimum is essentially impossible to find solely by a local optimization method. Therefore, a global optimization method should be used to solve this problem. The global optimization method used in this paper is a multilevel single-linkage (MLSL) method. It combines a stochastic method with a local optimization method, which has proved to be a very effective heuristic means of solving very complex global optimization problems. Therefore, axial permeability and transverse permeability of the core sample can be simultaneously determined with CDT method. We should note that the obtained permeabilities are apparent permeabilities if slippage effects exist. We mainly focus on the method for determining axial permeability and transverse permeability through CDT method. The absolute permeabilities can be obtained by a set of permeability tests under different gas pressure. Interested readers may refer to 26 for the details.

| Verification of the proposed method
To verify the accuracy of the proposed method, we performed numerical simulations of 3 cases that cover variations in permeabilities. The axial permeability and transverse permeability were calculated using the method described in Section 2.1, from the simulated gas production decay data, and then compare them to the values actually used in the numerical simulations.
In the numerical simulations, the core sample was discretized into 2500 grid cells (shown in Figure 3) and the canister was represented as a special grid cell attached to the sides of the core. The sample dimensions are representative of actual core tested in the laboratory flow tests, that is, 50 mm in diameter and 50 mm in length. These simulations were completed using TOUGH + REALGASBRINE (TOUGH+), a widely used Schematic diagram of numerical simulation model numerical simulator for nonisothermal multiphase flow (an aqueous phase and a real gas mixture) in a gas-bearing medium, with a particular focus in ultratight systems. 27 In the numerical tests, the temperature was set to 25℃, initial pressure and gas drawdown pressure were set to 3 MPa and 2MPa, respectively. Figure 4 shows the gas flow data generated from three numerical tests. With these data, we calculated the axial permeability and transverse permeability by history match. Note that there are small differences between the results from numerical simulations and analytical solution. This can be explained as follows. In the analytical solution, cumulative gas flow is calculated by the integral of gas flow from an infinitesimal area on the cylinder's periphery, but the simulator is not capable of modeling his process exactly. If the number of grids is infinity, numerical simulation should reproduce the analytical-solution results perfectly. This phenomenon is also observed in Hosseini's numerical tests, 24 which indicates that numerical-simulation results get closer to the analytical solution with the increase of grid number. Table 1 shows the comparison between the calculated permeability using the proposed method from the data generated by numerical simulations and the true values used in numerical simulation. As shown in Table 1, the determined axial permeability and transverse permeability values are very close to the true values. The relative error between the mean estimated value and the true value ranges from 1.7% to 4.5% for transverse permeability, and 1.9% to 5.3% for axial permeability, supporting the usefulness of our method. Figure 5 shows the gas production data versus time for the two samples. From the figure, the gas production increases sharply at the early stage due to the pressure gradient difference. The speed at which gas production increases slows down with flow time. Eventually, gas production stabilizes after the pressure attains the equilibrium pressure. The larger the initial pressure, the more the gas production. We then obtained the axial permeability and transverse permeability by history matching, as shown in Figure 6. Note that the accuracy of the fitting curves decreased with the increase of pressure. This effect arises from the adiabatic changes in temperature caused by changing pressure. Because all the CDT tests are conducted under room temperature and the gas drawdown pressure is set to 1 MPa, the increase of initial gas pressure will cause greater changes in temperature, which further influence the gas properties, leading to a larger fitting deviation. Figure 6 shows that all the permeability values parallel to bedding are bigger than the permeability values perpendicular to bedding. For example, the axial permeability of vertical sample is about 6 times of the transverse permeability, while this ratio decrease to about 0.2 for the horizontal sample. The permeability anisotropy calculated using our method lays very well within the permeability anisotropy values experimentally measured for the same formation and zone. 18

| Comparison between CDT results and PPD results
Pressure pulse decay (PPD) method is widely used for axial permeability measurement of tight reservoir rocks. To compare with CDT results, we performed three groups using the PPD method with the core samples. Figure 7 shows the schematic sketch of the PPD technique. [28][29][30][31] The sample is first in equilibrium and then opens the valve between the upstream chamber and the sample as well as the valve between the sample and the downstream chamber to let fluid flowing from the upstream chamber to the downstream chamber through the sample driven by the pressure difference. The axial permeability can then be estimated from the pressure pulse decay curve (the upstreamdownstream pressure difference with time) using mathematical solutions of the related flow problem.
Brace's solution is used to estimate axial permeability in PPD tests. 3 Table 2 shows the results of PPD tests on the tested samples. Note that transverse permeability cannot be obtained from the PPD method. Therefore, we only compare the axial permeability of the two methods. Compared with the permeabilities measured from CDT tests, the evolution of axial permeability versus gas pressure for the two methods has the same trend (shown in Figure 8). Specifically, permeability increases with the initial gas pressure. Axial permeability values obtained from the PPD method are less than permeability values of CDT method for the full range of gas pressure. The differences existing in the results of the two methods result from the following two reasons. First, PPD tests must be conducted under confining pressure to avoid lateral flow between the rock and the sleeve, which may close the microfractures in the rock and hence decrease the rock permeability. Second, though the initial gas pressures are almost the same for the two methods, the final equilibrium pressures are varied with each other, which causes different slippage effects on the measured results. Meanwhile, we estimated the permeability anisotropy using the axial permeability of the vertical and horizontal samples obtained from PPD tests. The permeability anisotropy values are 6.28, 7.45, and 7.68 respectively, slightly larger than the values obtained from CDT tests.

DISCUSSION
In this section, the effects of permeability anisotropy, porosity, initial gas pressure, rock dimensions, and adsorption on the profiles of gas production are investigated. The source of error and limitation of the proposed method are also explored. In the parametric study, we set up a base case and the parameters are as follows: k x = k r =1 × 10 -18 m 2 , ϕ = 5%, p i = 1.5MPa, L = 50mm, R = 25mm.

| Effects of permeability
The curves of cumulative gas volume versus time with variable permeabilities are shown in Figure 9 for cases, k x = k r , and Figure 10 for cases k x and k r are varied with each other. The figures clearly illustrate that gas volume produced in the canister is very sensitive to variation in transverse permeability and axial permeability. Both the transverse permeability and axial permeability affect the gas-flow transient time significantly, specifically, the lower the permeabilities are, the longer the gas-flow transient time is. However, these permeabilities have no effect on the cumulative gas volume.
For cases, k x and k r have different values, as shown in Figure 10, k r plays a more important role than k x in the evolution of cumulative gas volume versus time curve. For example, the gas-flow transient time of the cases for k r ranging from 1 × 10 -18 to 1 × 10 -19 m 2 increases much more than the cases for k x ranging from 1 × 10 -18 to 1 × 10 -19 m 2 . Due to that the cumulative gas volume is affected by both permeability and surface area of the core sample, the effects of k r and k x are also influenced by the aspect ratio of the core sample (ratio of the sample length to the diameter). When the aspect ratio is much less than 1 (ie, the sample length is less than the sample diameter), k x has more contributions to the cumulative gas volume than k r as a result of larger top-bottom surface area. With the same core sample volume, when the aspect ratio is

F I G U R E 5 Experimental results of CDT (A) for the vertical
sample. Experimental results of CDT (B) for the horizontal sample much larger than 1 (ie, the sample length is larger than the sample diameter), k r has more contributions to the cumulative gas volume than k x as a result of larger lateral surface area.

| Effects of porosity
The effects of porosity on the curve of cumulative gas produced versus time is shown in Figure 11. It is apparent that a large pore volume can store more gas and hence affect the gas-flow transient time. As shown in the figure, with the increase of porosity, gas production increases significantly, and the time required to attain final pressure equilibrium is also significantly extended. Figure 12 shows the effects of the initial gas pressure on the curves of cumulative gas volume versus time. As one would intuitively expect, a large initial gas pressure leads to more gas production but has no significant effect on gasflow transient time. While the initial gas pressure increases from 1.

| Effects of rock dimensions
To investigate the effects of rock dimensions, sample length and radius are investigated, and the results are shown in Figure 13. The figures indicate that larger rock size will increase both the cumulative gas produced and gas-flow transient time. In particular, the final cumulative gas produced shows a linear relation with sample length and quadratic function relation with sample radius, which indicates a linear relationship between the final cumulative gas produced and the sample volume. Therefore, it is advisable, in the interest of more gas production, to use a larger sample when testing a sample with small porosity.

| Effect of adsorption
Gas adsorption provides additional storage in the core sample. Therefore, desorbed gas from canister degassing tests should be taken into account if adsorbing gas (eg, CH 4 and CO 2 ) is used as testing gases for shale permeability measurement. Equation (1) can be revised as follow when accounting for desorbed gas 24 : where D is gas desorption source. To make the final partial differential equation analytically solvable, a linear adsorption relationship is adopted here. This assumption is reasonable at low pressures. 11 where ψ is the amount of gas desorption, and γ is linear adsorption parameter, which can be estimated from actual Langmuir adsorption parameters.
where V m and p L are Langmuir volume and pressure constants. Based on the linear adsorption assumption, D can be calculated as follows: where ρ s is rock density. Substituting Equation (19) into Equation (16), and replacing the dimensionless parameters, we obtain: where = s c The boundary and initial conditions are the same to Equation (7). With the same analytical solving method to Section 2.2, we obtain the total gas production coupled with gas desorption as follows: Figure 14 shows the gas production curves with and without consideration of gas desorption (p L = 2MPa, V m = 2.8 cm 3 /g). As expected, the amount of gas production with gas desorption is much larger than the case without gas desorption. The difference is more significant for rocks with larger adsorption capacity. Meanwhile, due to the effect of gas desorption, the gas-flow transient time also increases. We should note that low gas pressure should be used when testing adsorptive rocks due to the assumption of linear adsorption. In addition, gas desorption may induce matrix shrinkage which potentially increases the rock permeability. 32,33 For these reasons, it is advisable to use nonadsorptive gases when testing rock permeability with CDT method.

| Source of error and limitation
Several sources of error exist in the CDT method, such as system leakage, the measurement device, temperature fluctuation, errors from the determination of pore volume, and so on. Temperature fluctuation arises from the adiabatic changes in temperature caused by suddenly changing in gas drawdown pressure. To eliminate such effect, a small pressure difference between the saturation pressure and the gas drawdown pressure should be applied. Since the CDT method records the gas flow data, a flowmeter with high resolution is needed in case of testing permeabilities of rocks with small pore volume and ultralow permeabilities. Moreover, the proposed method of estimating the axial permeability and transverse permeability is an inverse problem which may have nonunique solutions; therefore, reasonable initial values of axial permeability and transverse permeability should be given to match the produced gas data (eg, the permeability (19) D = − s dp dt