Analytical models for determining pressure change in an overlying aquifer due to leakage

Available online at www.sciencedirect.com

Energy Procedia 4 (2011) 3833–3840

Energy Procedia www.elsevier.com/locate/procedia

Energy Procedia 00 (2010) 000–000 www.elsevier.com/locate/XXX

GHGT-10

Analytical Models for Determining Pressure Change in an Overlying Aquifer due to Leakage Mehdi Zeidouni1*, Mehran Pooladi-Darvish, David W. Keith Department of Chemical and Petroleum Engineering,University of Calgary, 2500 University Drive, Calgary T2N 1N4, Canada Elsevier use only: Received date here; revised date here; accepted date here

Abstract Various methodologies are proposed to reduce CO2 emissions that are believed to be the main drivers of the climate change. CO2 capture and storage in deep underground formations is one of the promising methods that allow reducing the emissions while continuing the use of fossil fuels. Injection of immense quantities of CO2 is required to make a reasonable cut of the emissions. Deep saline aquifers can provide the capacity to accommodate the storage of such huge amounts of CO2. However, one of main challenges in deployment of CO2 storage is the risk of CO2 leakage through pathways in the cap-rock overlying the target aquifer. The sealing capacity of the cap-rock must be evaluated to ensure the safety of the storage. Therefore, characterization of the cap-rock is required to find the potential leakage pathways even before the CO2 storage begins. Methods to characterize the leakage pathways are proposed at two different scales: 1- by point sampling of the cap-rock and testing the potential pathways such as abandoned wells,2- by analysing geophysical (e.g. 3-D seismic) data to estimate paths of upward migration of the injected CO2. Flow based methods have the potential for bridging the large gap that exists between the length scale of these two approaches. The aquifer could be tested for the leakage pathways before CO2 storage. This will allow finding proper storage aquifers and locations for the injection wells. In this work we present an analytical model to evaluate the pressure variation in the overlying aquifers due to leakage from the storage aquifer. In a companion paper, this model will be used along with an inverse modelling approach to locate and characterize the leakage pathways based on pressure data. This paper introduces two new analytical solutions: 1- exact solutions for the pressure variation in an overlying aquifer due to leakage (obtained in Laplace-transformation domain), 2- time-domain approximations for the exact solutions to make the inversion possible. In deriving the analytical solutions two aquifers are considered: storage and monitoring. The aquifers are separated by an aquitard and are in communication through a leakage pathway. In departure from previous works the leakage pathways are not required to be line source/sink. Such consideration allows incorporation of large pathways such as stratigraphic and structural heterogeneities in the cap-rock. We consider a single-phase 1-D radial flow system in the storage and monitoring aquifers. Both of the aquifers are considered as homogeneous, isotropic, and infinite-acting with constant thickness. The injection (or production) rate is taken as constant. The analytical solution are applied to a base case and corroborated versus numerical solution. ©c 2010 Ltd.byAll rights reserved ⃝ 2011Elsevier Published Elsevier Ltd. "Keywords: Pressure monitoring; Analytical solutions; Leakage"

* Corresponding author. Tel.: +1-403-220-7929; fax: +1-403-210-3894. E-mail address: [email protected].

doi:10.1016/j.egypro.2011.02.319

3834

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

1. Introduction Capture of CO2 from industrial point sources and long-term storage in underground geological formations provides a feasible option to cut the CO2 emissions while maintaining access to fossil energy. The number of projects involving the underground storage of CO2 is increasingly growing. The current scale of the projects ranges from industrial operations -around 1 million tons (Mt) of CO2 injection per year- to much smaller-scale research projects – injecting thousands of tons per year- the goal of which is to learn more about reservoir dynamics associated with CO2 injection. In order to cut the emissions noticeably the size of the projects must increase to much larger scales as Gt-scale. Deep saline aquifers are believed to be able to provide the required capacity for such immense quantities of CO2. Unlike oil and gas reservoirs where holding hydrocarbons for millions of years provide an evidence of good isolation capacity by an impermeable cap-rock, such observations are missing for saline aquifers. The aquifers may be overlain by a leaky seal which allow migration of the CO2 to the upper aquifers and may lead to leakage to the surface. To accommodate safe storage of CO2 characterization of the cap-rock is required by locating and characterizing the leakage pathways. The operator may then decide to seal the pathways or inject the CO2 far away from the leaks. The potential leakage pathways could include: a) transmissive faults and fractures, b) abandoned wells (penetrating entire seal or part of it), c) active wells that partially penetrate the seal, d) and local seal weaknesses. Local weakness in the cap-rock can be a channel dug through the shale by gravitationally unstable deposits. This leaves coarse-grained, high permeable channels and lobes in the cap-rock through which CO2 can leak to shallower formations [1]. Recent time-lapse seismic data of the Sleipner project have shown that there may be such a sand-fill in the intra-reservoir mudstone. As a result of broken continuity the mudstone did not significantly impede the CO2 vertical flow in Utsira formation. The channel may be of meter-scale diameter based on the size of leaked chimney [2]. Also, karstic collapse due to natural voids or weaknesses in the rock can initiate cap-rock voids [3]. In carbonate environments such as Wabamun Lake Sequestration Project (WASP) karsting features as large as 0.2-1 km2 may exist in Calmar cap-rock [4]. Migration through the cap-rock causes a pressure change in the upper aquifers. Such pressure change is of interest because it can be used to infer the leakage pathway location and characteristics. The analytical model to be developed here is not appropriate for leakage through faults. The goal of this study is to develop analytical models suitable for leakage from abandoned wells and/or other local weaknesses in the cap-rock. In this work we present an analytical model to obtain the pressure change caused by leakage from the storage aquifer to an upper aquifer. The pressure change is first obtained in Laplace-transform domain. Using the late-time approximation of the functions involved in the solution an asymptotic solution is obtained. Equations to evaluate the pressure change at the location of monitoring well and the leakage rate are also derived. The exact and asymptotic analytical solutions are then compared by applying to a base case. The solutions are also corroborated with the numerical solution provided by commercial simulator. Study of potential leakage of hazardous wastes into shallow groundwater aquifers has started more than two decades ago [5]. Several analytical and semi-analytical models have been developed on leakage detection and quantification. Javandel et al. (1988) developed an analytical model for leakage from a horizontal aquifer to an upper aquifer through an abandoned well. The upper aquifer is separated from the lower aquifer -in which water is injected at a constant rate- by an impermeable aquitard. The pressure at the storage aquifer is monitored through a monitoring well. A method to calculate the amount of leakage from an abandoned well and the corresponding drawdown at the monitoring well(s) is obtained. The model considers single-phase flow with constant density. The wells (injection and leaky) are assumed as line source/sink that fully penetrate the aquifers. The pressure is assumed to be constant in the upper aquifer. Both the pressure drawdown at the monitoring well due to leakage and the leakage rate are described by analytical expressions which involve integrals to be numerically evaluated. Avci (1994) have solved this problem by elimination of the assumption of constant pressure at the upper aquifer [6]. As a result lower leakage rate was obtained. He also solved the above problem for leakage through an improperly plugged well which allow communication between the aquifers. The obtained solutions are in Laplace transform domain. The same problem under the steady-state conditions has been solved earlier by Silliman and Higgins (1990) for fully penetrating wells and by Avci (1992) for partially penetrating wells [7, 8]. Considering the behaviour of the

3835

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

leakage flux (continuously increases, rapidly at the beginning and then at a continually decreasing rate) Nordbotten et al. (2004) approximated the leakage rate[9]. They obtained a real-time approximation for the leakage rate. In departure from Avci’s work [6] the exact analytical solution in this work considers the leak to have a finite size (not line source/sink). In future works, the analytical solutions obtained in this work will be used for leakage detection and characterization. Note that the geomechanical effects are not included in our analysis. Such effects can cause uplifting of the cap-rock and even fracturing. This may lead to pressure change in the upper aquifer also. 2. Physics of leakage problem The physical configuration and the variables used in the problem formulation are shown in Figure 1. Consider two aquifers: target storage aquifer to be tested for appropriateness for CO2 storage, and upper monitoring aquifer where the pressure is to be monitored. A single-phase 1-D radial flow system is considered in the aquifers which are separated by an impermeable aquitard. The leakage occurs in vertical direction through a single leakage pathway of radius rl and permeability kl. The leakage pathway is at a distance R from the injection well and distance U from the monitoring well. Knowing that the aquifers are very large the time to the end of the infinite-acting flow mechanism is assumed to be very large. However, if the aquifer is limited by a boundary e.g. a sealing fault one can use the image method to count for the effect of such boundary on the leakage rate using the infinite-acting solution. The aquifers are also considered isotropic and homogeneous with known and constant properties (e.g. permeability, porosity, thickness, compressibility). The injection fluid is injected at a constant rate q and considered to have same properties as the aquifer brine. The aquifers are not in communication with a third aquifer.

3. Analytical model To obtain the pressure variation and leakage rate we decompose the system into 4 problems (components) and then combine the results. The starting point is the pressure variation in the monitoring aquifer in response to unknown and time-dependent leakage rate (ql) which is the solution to the pressure diffusivity equation centered at the leakage pathway. Based on Darcy’s equation, the leakage rate is a function of the pressure difference between the two aquifers at the location of the leak. To obtain the pressure at the location of the leak in the storage aquifer the superposition principle can be used. The pressure response to injection can be obtained by solving the diffusivity equation centered at the injection well under constant rate boundary condition. Pressure response to leakage is achievable through solving the diffusivity equation considering the leakage path as center. Superposition of the two solutions evaluated at the location of the monitoring well provides an equation for the pressure variation at the monitoring well location in the storage aquifer. Combination of such equation with pressure variation in the monitoring aquifer gives an equation for the time-varying leakage rate. Details of deriving the problem solution can be found in [10]. The following is obtained for the pressure change at the monitoring well:

§

qlD ( s ) K 0 ¨ PmD TD rlD

Where:

s

·

UD ¸

© KD ¹ § s · s K1 ¨ rlD ¸ KD © KD ¹

Equation 1

3836

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

K0 3 2

s



K  s s RD

1

qlD

§

K 0 ( srlD ) rlD

s K1 ( srlD )

·

s

K0 ¨

© KD

Equation 2

rlD ¸

¹  § s · D K1 ¨ rlD ¸ KD © KD ¹



1

s

TD rlD

and:

R

RD

rw

2S k s hs

PmD

U

UD

,

qP B

P

mi

rw

L

LD

,

Kst

 Pm , t D

2

rw

qlD

,

rl

rlD

,

,

rw

km hm

TD

ks hs

,

KD

Km Ks

2

,D

rl kl

,

2ks hs hl

ql

q Attempts for analytical inversion of the above solutions into time domain was not successful. These expressions can be numerically inverted using the method developed by Stehfest [11]. In the following, we present asymptotic solutions in real-time domain. rw

4. Asymptotic solution

A late-time asymptotic solution for pressure change at the monitoring well and the leakage rate are obtained the detailed derivation of which is given in [10]. Using Bessel functions’ properties we simplify Equation 1 and Equation 2. We then use the approximate inversion of some of the terms that appear in the simplified solution to get the following real-time approximate solutions for pressure change at the monitoring well and the leakage rate respectively:

§ J 1 § RD2 · · qlD § § U2 ·· N  ln ¨ D ¸ ¸   ln ¨  ¨ ¸ ¸ ¨ 1  TD © 2 2 © 4t D ¹ ¹ 2TD © © KD ¹ ¹

Equation 3

1  C  t N  ln  R

Equation 4

1

PmD

1

qlD

2

1  1 / T

D

D

D

Where: J 

C  tD

1

J



N  2J  ln  4t N  2J  ln  4t D

D

§  T 1 ¨K  ln D 2 ¨ r D  TD  1 ¨ lD © 1

N

2TD

D

2



S

2

J 

S

2

J  2[ (3) 6 2  3 N  2J  ln  4tD  N  2J  ln  4tD 4 2

3

· ¸ ¸ ¸ ¹

Equation 5

Equation 6

and[(3) =1.2020569032 is the Riemann Zeta function. Equation 3 gives the dimensionless pressure at the monitoring well, where measurement will be made. 5. Numerical validation and comparing the exact and asymptotic solutions

In this section the exact and asymptotic analytical models are compared while validated through numerical investigations. The leakage rates, the pressure change at the monitoring well, and that in the storage aquifer are

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840

3837

M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

compared using the analytical and numerical solutions. The numerical solutions are obtained using commercial Black-Oil reservoir simulator [12] which can easily handle a single-phase flow. The numerical and analytical solutions are applied to a base case the properties of which are given in Table 1. For the base case water in injected at constant rate of 864 m3/day in a 30 m thick storage aquifer. The upper aquifer is 10 m thick which is separated from the storage aquifer by a 4 m thick cap-rock. The permeability of the storage and monitoring aquifers are 100 and 25 mD respectively. A leakage pathway is located at 50 m distance from the injection well and 35 m distance from the monitoring well. The radius and permeability of the leak are 0.2 m and 150 mD respectively. The base case is modelled numerically using IMEX black-oil simulator to check the numerical versus the analytical results. The cylindrical geometry of the leaky path cannot be accurately modelled using the Cartesian grids. However, HYBRID grid option in IMEX allows refining a Cartesian grid block into a cylindrical grid in the centre plus an outer block whose outer boundary fits the shape of the parent grid. We use this option to properly model the leaky pathway. In the cap-rock formation, the outer grid is NULLed by assigning a porosity of zero. An area of 1.135*1.135 m2 is considered for the block where the leaky well is located (since the leak diameter is automatically taken to be 0.35 of the grid block side). The exact analytical solution is compared to the asymptotic and numerical solutions in figures 2 and 3 in terms of the monitoring pressure and leakage rate. The monitoring duration is considered to be 100 days. The leakage rate over the first three months is approximately 4 orders of magnitude smaller than the injection rate. The pressure change in the monitoring well in the monitoring aquifer is less than 3 kPa. A very good match between the exact, asymptotic, and the numerical solutions is observed. 6. Conclusions

An analytical model is presented to estimate the leakage rate and pressure change due to leakage at the storage and monitoring aquifers. The model is obtained by solving the flow equations in the monitoring and storage aquifers which are coupled by the flow rate at the leak. Superposition principle is used to obtain the pressure change at the location of the leak in the storage aquifer. The exact analytical solution is obtained in Laplace domain. Real-time asymptotic solution is obtained through approximation of the functions involved at late-time period. The solutions are corroborated against the numerical simulation results for a base case. Close agreement between the analytical and numerical results is found. In future work, the analytical solution developed here will be used to infer the leak location and properties based on pressure data measured at a monitoring well.

Acknowledgements

Research leading to this note has been supported by a Natural Sciences and Engineering Research Council of Canada (NSERC) Strategic Grant and by the Alberta Energy Research Institute (AERI), with additional funding from industry partners. The first author thanks Chris Eisenger for discussions on subsurface geology. Nomenclature: A B h K0 K1 k L P q R r s

cross sectional area, m2 water formation volume factor, vol. @ Res. cond./ vol. @ St. cond. thickness, m Zero order modified Bessel function of the second kind First order modified Bessel function of the second kind permeability, m2 Monitor-injector distance Pressure, Pa volumetric flow rate, m3/s leak-injector distance, m radius, m Laplace transform dummy variable

3838

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

[(3) 

Transmissivity=permeability × thickness, m3 time, s leakage coefficient Difference from the initial value Euler constant = 0.5772… aquifer diffusivity coefficient Permeability / (Porosity × fluid viscosity × total compressibility) , m2/s viscosity, Pa.s leak-monitor distance 1.2020569032, Riemann Zeta function

Subscripts: D i l m s w

dimensionless initial leakage monitoring aquifer storage aquifer well

Superscripts: prime: ‘

shows the radius originated at the injection well

T t

D ' J K P U

References: 1. Grimstad, A.-A., Georgescu, S., Lindeberg, E., Vuillaume, J.-F., 2009, Modelling and Simulation of Mechanisms for Leakage of CO2 from Geological Storage, Energy Procedia 1, 2511-2518. 2. Chadwick R. A., Arts R., Eiken O., 4D seismic quantification of a growing CO2 plume at Sleipner, North Sea, Geological Society, London, Petroleum Geology Conference series 2005, v. 6, p. 1385-1399, doi: 10.1144/0061385. 3. Waltham, T., Bell. F. G., Culshaw, M. G., 2005, “Sinkholes and Subsidence: Karst and Cavernous Rocks in Engineering and Construction”, Springer Berlin Heidelberg. 4. Alshuhail A., Lawton D., Isaac H., 2009, Wabamun Area CO2 Sequestration Project (WASP): Seismic Characterizations of the Nisku Formation, available from: http://www.ucalgary.ca/wasp/reports.html 5. Javandel, I., Tsang, C.F., Witherspoon, P.A., Morganwalp, D., 1988, Hydrologic Detection of Abandoned Wells Near Proposed Injection Wells for Hazardous Waste Disposal. Water Resour. Res. 24. 6. Avci, C. B. , 1994, Evaluation of flow leakage through abandoned wells and boreholes, Water Resour. Res., 30(9), 2565–2578. 7. Silliman, S., Higgins, D., 1990, An Analytical Solution for Steady-State Flow Between Aquifers Through an Open Well. Ground Water 28, 184-190. 8. Avci, C.B., 1992, Flow occurrence between confined aquifers through improperly plugged boreholes. Journal of Hydrology 139, 97-114. 9. Nordbotten J. M., M. A. Celia, S. Bachu, 2004, Analytical solutions for leakage rates through abandoned wells, Water Resources Research, 40(4), W04204, doi:10.1029/2003WR002997. 10. Zeidouni M., Pooladi-Darvish M., Keith D., 2010. Inter-Aquifers Flow through Leakage Conduits: Analytical Solutions for Leakage Rate and Pressure Change Evaluation, Submitted 11. Stehfest, H ., Algorithm368 numerical inversion of Laplace transforms, D-5, Commun. ACM, 13(1), 47-49, 1970. 12. IMEX Version 2008 user’s guide, Computer Modelling Group Ltd, Calgary, Alberta, Canada, http://www.cmgl.ca/software/imex.htm

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000 Table 1. Descriptions of the base case and lower and upper limits for sensitivity analysis

kl, leak permeability (m2)

150e-15

R, leak-injector distance (m)

50

U, leak-monitoring distance (m)

35

hl, leakage interval (m)

24

rl, leak radius(m)

0.2

P, brine viscosity (Pa.s)

0.5e-3

monitoring aquifer compressibility (1/Pa)

1e-9

monitoring aquifer porosity(fraction)

0.1 2

km, monitoring aquifer permeability(m )

25e-15

hm, monitoring aquifer thickness(m)

10

Storage aquifer compressibility (1/Pa)

1e-9

Storage aquifer porosity (fraction)

0.1 2

ks, Storage aquifer permeability (m )

100e-15

hs, Storage aquifer thickness(m)

30

q, injection rate into the storage aquifer (m3/s)

0.01

rw, injection well radius(m)

0.1

B, water formation volume factor

1

Figure 1. Schematic view of the aquifer-leak-aquifer system

3839

3840

M. Zeidouni et al. / Energy Procedia 4 (2011) 3833–3840 M. Zeidouni et. al./ Energy Procedia 00 (2010) 000–000

3000 2500

'Pm (Pa)

2000 Exact Asymptotic Numerical

1500 1000 500 0 0

20

40

Time(day)

60

80

100

Figure 2. Pressure change at the monitoring well versus time for the base case

2.5E06

Leakagerate(m3/s)

2.0E06 1.5E06

Exact Asymptotic Numerical

1.0E06 5.0E07 0.0E+00 0

20

Figure 3. Leakage rate versus time for the base case

40

Time(day)

60

80

100