Geomechanical stability of CO2 containment at the South West Hub Western Australia: A coupled geomechanical–fluid flow modelling approach

International Journal of Greenhouse Gas Control 37 (2015) 12–23

Contents lists available at ScienceDirect

International Journal of Greenhouse Gas Control journal homepage: www.elsevier.com/locate/ijggc

Geomechanical stability of CO2 containment at the South West Hub Western Australia: A coupled geomechanical–fluid flow modelling approach Y. Zhang a,∗ , L. Langhi b , P.M. Schaubs a , C. Delle Piane b , D.N. Dewhurst b , L. Stalker b , K. Michael b a b

CSIRO Mineral Resources National Research Flagship, PO Box 1130, Bentley, WA 6102, Australia CSIRO Energy, PO Box 1130, Bentley, WA 6102, Australia

a r t i c l e

i n f o

Article history: Received 14 November 2014 Received in revised form 7 February 2015 Accepted 4 March 2015 Keywords: South West Hub Geomechanical stability CO2 injection and storage Geomechanical and fluid flow modelling Lesueur Formation

a b s t r a c t An area in the Southern Perth Basin has been selected as a potentially suitable site for CO2 injection and storage as a part of the South West Hub Project (SW Hub), due to its proximity to major CO2 emission sources and the presence of potentially suitable geology. This 3D modelling study attempts to assess the geomechanical stability of faults and intact host rocks during CO2 injection at the SW Hub. The stratigraphy and fault structure of the 3D model are based on the architecture of an E–W cross section in a pre-existing 3D geological model that represents a comprehensive synthesis of seismic, stratigraphic and structural data. In the models, the rocks and faults are simulated as Mohr–Coulomb elastic–plastic materials, and their geomechanical and hydrological properties are based on experimental data from the Harvey-1 drill core samples and also information from literature. A series of models are performed to assess five injection scenarios with injection rates of 1–5 million tonnes per year over a period of 20 years. The results show that the simulated CO2 injection scenarios would not lead to fault reactivation or breach the overlying Yalgorup or Eneabba Shale formations in the area. Some small smooth uplifts are recorded as a result of injection. In the models assuming weak faults, average ground surface uplifts are 0.4–1.8 cm for the injection rates of 1–5 million tonnes per year, over an area of approximately 2.5 km radius around the hypothetical injection site. Uplifts are marginally smaller when assuming strong faults. Crown Copyright © 2015 Published by Elsevier Ltd. All rights reserved.

1. Introduction An area in the southern part of the Perth Basin in Western Australia was originally proposed by the Australian Petroleum Cooperative Research Centre’s GEODISC project as an environmentally suitable site for CO2 injection, capture and storage (e.g. Bradshaw and Rigg, 2001), due to its proximity to major CO2 emission sources and the presence of potentially suitable geology (see also Causebrook et al., 2006; Varma et al., 2009). A decision was subsequently made to proceed with a comprehensive site assessment and geological modelling of the sub-basin structural unit known as the Mandurah Terrace. This site is referred to as the South West CO2 Geosequestration Hub (SW Hub; Fig. 1). The effort of new data acquisition in the area has included a 2D seismic survey in 2011, the

∗ Corresponding author. Tel.: +61 8 6436 8626; fax: +61 8 6436 8555. E-mail address: [email protected] (Y. Zhang). http://dx.doi.org/10.1016/j.ijggc.2015.03.003 1750-5836/Crown Copyright © 2015 Published by Elsevier Ltd. All rights reserved.

drilling of the Harvey-1 stratigraphic well (total depth was 2945 m MDRT [measured depth from rotary table]) in 2012 and associated coring and logging of well data, as well as subsequent laboratory geomechanical tests and structural analyses. A comprehensive discussion and summary of previous key studies for the SW Hub and the evolution of the geological understanding of the Lesueur sandstone and the other associated formations in the stratigraphic sequence of the area can be found in Stalker et al. (2013) (also see Crostella and Backhouse, 2000; Bradshaw et al., 2000; Hennig and Otto, 2005; Ennis-King and Wu, 2005; Causebrook et al., 2006; Dance and Tyson, 2006; Varma et al., 2007, 2009). A major conclusion from these previous studies is that the SW Hub does have the potential for being a commercial scale CO2 capture and storage (CCS) site, subject to confirmation of the CO2 trapping mechanisms. Delle Piane et al. (2013) carried out a comprehensive study characterizing the faciesbased rock properties and storage capacity of the geological units intercepted along the Harvey-1 well. In particular, their laboratory

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

13

Fig. 1. (a) Regional geological map showing the location and structural outline of the South Perth Basin and location of the South West Hub (SW Hub). (b) The structural map of the SW Hub area showing fault traces and dip directions for the area in the dashed line box in (a). The A and B line defines the trace of an E–W cross section in the SW Hub area. Both maps are after Langhi et al. (2013).

geomechanical tests generated a new data set of density, elastic parameters, cohesion, friction angles, porosity and permeability values for the Yalgorup member and Wonnerup member of the Lesueur sandstone, the two key rock units of the Harvey-1 well stratigraphic sequence. The geological and structural framework of the SW Hub area was constructed by Langhi et al. (2013). Their initial model is based on the interpretation of 2D seismic reflection data for the geometries of key faults and stratigraphic units and also integrates stress data (Pevzner et al., 2013), facies and rock property data (Delle Piane et al., 2013; Olierook et al., 2014) and a 3D facies model (Griffiths et al., 2012). This multi-disciplinary analytical effort led to the construction of a 3D geological model for the area, characterizing the fault and stratigraphic network of the region (Figs. 1 and 2). The geological data and knowledge for the region described above form the necessary basis for advancing our understanding of the carbon storage potential at the SW Hub by gaining an insight into its geomechanical behaviour under CO2 injection conditions. Several previous numerical modelling studies for other CCS sites in the world (e.g. Rutqvist et al., 2010; Vidal-Gilbert et al., 2010; Verdon et al., 2013; Rinaldi et al., 2013, 2014; Aruffo et al., 2014; Mbia et al., 2014) have demonstrated the successful applicability of a numerical modelling approach to the assessment of the geomechanical and fluid flow behaviours of such sites. In this paper, we present the results of a group of coupled geomechanical–fluid flow models with one-phase fluid (simulated as water)

attempting to assess the geomechanical stability of the faults and host rocks during CO2 injection at the SW Hub. In addition, we examine the potential deformation (ground surface uplift) at the site for a number of injection rates and rock/fault geomechanical property scenarios. 2. Geological background The SW Hub area (about 110 km south–southwest of Perth, Western Australia) is located in the central part of the onshore Perth Basin at the southern end of the Mandurah Terrace, east of the off-shore Vlaming Sub-basin (Fig. 1a). At a regional scale, the area is bounded by the N–S trending Darling fault to the east and by the Badaminna fault and Harvey Ridge to the west and southwest, respectively. The stratigraphy of the Central Perth Basin is composed of sedimentary rocks with ages from the Permian to Cretaceous (e.g. Crostella and Backhouse, 2000; Stalker et al., 2013). Fig. 2 gives the details of the stratigraphic units intercepted at the Harvey-1 well. The most important stratigraphic units in the sequence are described here. The Wonnerup member of the Lesueur sandstone is a predominantly coarse to gravelly sandstone unit, overlying the relatively thin Sabina sandstone. This unit represents the key aquifer unit in the area suitable for CO2 storage. The Yalgorup member is mainly an inter-bedded package of mixed coarse to gravelly sandstone, fine to medium grained sandstone, siltstones and

14

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

Fig. 2. Stratigraphic section for rock units intercepted at the Harvey-1 well (after Langhi et al., 2013). FM–Formation.

mudstone. The Basal Eneabba shale unit and the Yalgorup unit are regarded as the containment formations at the SW Hub. The Eneabba Formation overlying the basal shale consists of fluvial coarse sandstone inter-bedded with local minor conglomerate, claystone and siltstone (Mory and Iasky, 1996). The structural analysis of Langhi et al. (2013) shows that the fault population in the SW Hub area is characterised by a series of mainly WNW-, NW- and NNW-trending small faults (F2 , F4 , F5 , F6 , F11 –F15 ) bounded by two larger faults (the N–S to NNWtrending F1 and NW to NNW trending F10 ) (Fig. 1b). The Harvey-1 well is located between the F5 and F9 faults, both ENE dipping. The stress field of the region is dominated by approximately E–W oriented maximum principal stress ( hmax ) and N–S oriented minimum principal stress ( hmin ), with the intermediate principal stress ( v ) being vertical (Hillis and Reynolds, 2000; Rasouli et al., 2013; Pevzner et al., 2013). This stress pattern is adopted as the initial stresses of our models. More specifically, a gravitational stress gradient is specified for the vertical stresses in our models, which are defined as  v = density × gravity × depth. The corresponding vertical stress gradients vary between 22.2 and

23.5 MPa/km due to structural and sediment density variations. The horizontal stresses ( hmax and  hmin ) are specified according to the ratios of  hmax / v = ∼1.203 and  hmin / v = 0.829 ( hmax gradient = 26.7–28.3 MPa/km and  hmin gradient = 18.4–19.5 MPa/km; see Langhi et al., 2013). All the stresses are compressive. Such stress ratios will lead to relatively small  hmax and small differential stresses and Mohr circles for shallow reservoir levels as will be illustrated subsequently. This is consistent with the magnitudes of  hmax and differential stress reported for the southern Perth Basin (e.g. Van Ruth, 2006; Pevzner et al., 2013). 3. Numerical model 3.1. Brief description of numerical method The present numerical modelling study has used a 3D finite difference code, FLAC3D (Fast Lagrangian Analysis of Continua; Cundall and Board, 1988; Itasca, 2005), which is capable of simulating the interactions between deformation and fluid flow in porous media. A 3D mesh consisting of hexahedral elements and

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

representing rock materials has been generated to fit the geometries of the geological structures to be simulated. Each element in the mesh behaves according to prescribed geomechanical and hydraulic laws in response to the applied boundary conditions. The material can yield and deform plastically and the mesh deforms and moves with the material. The code has been applied to the fields of structural geology (e.g. Ord, 1991; Strayer et al., 2001; Zhang et al., 2008), economic geology (e.g. Sorjonen-Ward et al., 2002; McLellan et al., 2004), petroleum geology (e.g. Langhi et al., 2010; Zhang et al., 2009, 2011) and CO2 storage simulation (e.g. Rutqvist et al., 2010; Röhmann et al., 2013). For the present study, rocks are simulated as Mohr–Coulomb elastic–plastic materials (Cundall and Board 1988; Itasca, 2005; Ord, 1991). The necessary constitutive parameters include shear modulus, bulk modulus, cohesion, tensile strength, friction angle and dilation angle. Under mechanical loading, the material deforms initially in an elastic manner up to a yield point and then deforms plastically upon yield. The yield is governed by the Mohr–Coulomb yield criteria, that is, the occurrence of yield when the maximum shear stress reaches a threshold magnitude. In addition, tensile failure occurs when the effective minimum principal stress is in tension and overcomes rock tensile strength. Dilation (positive volume change) can occur in the model with shear deformation. The dilatant potential of the Mohr–Coulomb material for plastic deformation is characterised by the dilation angle. FLAC3D adopts the engineering stress-sign convention, that is, compressive stress is assumed negative and tensile stress is positive. Single phase fluid flow is governed by Darcy’s law (Itasca, 2005; also see Mandl, 1988) for an isotropic porous medium and is coupled with geomechanical deformation in the model. Hydrogeological parameters including fluid density, fluid bulk modulus, porosity and permeability are assigned and pore fluid pressure is initialised before the simulation commences. Fluid flow velocities

15

are primarily a function of gradients in pore fluid pressures and variations in permeability. The interactions between mechanical deformation and fluid flow are reflected in three ways: (1) volumetric strain or solid volume changes result in pore pressure changes (for example, a local volume increase leads to a local pore pressure decrease); (2) changes in pore pressure cause changes in effective stress, which further affect the geomechanical response of the solid (e.g. a reduction in effective stress may induce plastic yield); and (3) the formation of any topographic elevation or depression in a model, as a result of bulk deformation, can change pore pressure or head distribution, and then lead to changes in fluid flow patterns. 3.2. Model geometry and property The stratigraphy and fault structure of the present model are based on the architecture of an E–W cross section (Fig. 4a) in the 3D geological model of Langhi et al. (2013) (Fig. 3), which trends through the Harvey-1 well location (section A and B in Fig. 1b). The model has been constructed in 3D, however, to allow out of plane fluid flow and ground movements by the projection of the 2D sectional geometries in the N–S direction (Fig. 4b). The resulting 3D model is 20.08 km in the E–W direction, 12.1 km in the N–S direction and 5 km in depth. Four faults (F1 , F9 , F5 and F10 from east to west) are included in the model, perpendicular to the E–W cross section, and are simulated as narrow zones in continuum contact with host rocks (i.e. fault contact planes with surrounding rock media are not simulated as discrete contacts). Fault dip angles vary between 52◦ and 75◦ . In the model, the stratigraphy of the SW Hub area is simplified into 7 rock units (Table 1; also see Fig. 2). The rock unit at the top of the model (Fig. 4) combines the Guildford and Leederville formations for simplicity, as both are thin stratigraphic units and

Fig. 3. The architecture of the 3D geological model of Langhi et al. (2013). Only four stratigraphic units are shown here, including, from top to bottom, Basal Eneabba shale, Yalgorup member, Wonnerup member and Sabina sandstone (see Fig. 2). The unit for dimension is metres.

16

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

Fig. 4. (a) Fault and stratigraphic geometries on an E–W cross section in the SW Hub area (A and B in Fig. 1b). The lower end of the dashed line shows the CO2 injection depth level at approximately the Harvey-1 well location. (b) The geometries and sizes of the 3D numerical model based on the 2D cross section above.

have comparable lithology (sandy clay with local conglomerates or sandstone with shale inter-beds and lignite seams). Below the Guildford–Leederville cover, all the other stratigraphic units dip east and form half-graben structures at the F1 and F10 faults. As stated above, the Wonnerup unit is the main aquifer unit, and lies underneath the major containment formations of the Yalgorup and Basal Eneabba shale units. The mechanical and fluid flow properties of the model are listed in Table 1. The properties for the Wonnerup and Yalgorup units are based on the average values from the laboratory experimental work on the drill core samples of the Harvey-1 hole (Delle Piane et al., 2013; Olierook et al., 2014); core samples are only available for these two units at Harvey-1. The parameters for the Basal Eneabba shale are based on the shale data used in Langhi et al. (2010) as core recovery in this interval was poor and no mechanical tests were possible. The values were chosen as shales usually display strengths lower than sandstones and these values are typical of shale strengths measured from other basins (e.g. Dewhurst and

Hennig, 2003). The properties of other units are chosen based on the dominant rock types (i.e. predominantly sandstone or sandstone inter-bedded with siltstone, shale and claystone) and information from the literature (e.g. Sarda et al., 1993; Plumb, 1994; Chang et al., 2006). The specification of tensile strengths in the model (i.e. being half of cohesion for each rock unit) is a simplified approach due to lack of data. This is consistent with general expectation that rock tensile strengths are lower than cohesions (e.g. Price and Cosgrove, 1990; Fang and Daniels, 2006), and in particular, with the results of Vanicek’s triaxial and tension tests (2013), which showed the tensile strength is about half of the cohesion of the tests. In the absence of direct laboratory measurements, a dilation angle of 2◦ is assumed for all the rock units. Dilation angle is the parameter determining volume increase (dilatancy) in the modelled Mohr–Coulomb elastic–plastic material when plastic shear strain is involved (e.g. Ord, 1991). This parameter is unimportant in the present model because there is no plastic failure under the present model conditions as will be shown below.

Table 1 Material properties of the model. Unit

Density (Kg m−3 )

Young’s modulus Poisson’s (GPa) ratio

Leederville Fm and Guildford Fm Eneabba Fm Basal Eneabba shale Yalgorup member Wonnerup member Sabina sandstone Willespie Fm Weak fault

2265

14

0.26

2265 2400 2260 2269 2269 2265 2280

14 2 7.9 20.3 20.3 14 2

0.26 0.35 0.23 0.295 0.295 0.26 0.15

7 2.75 5.28 14.1 14.1 7 2

Strong fault

2280

40

0.3

28

Cohesion (MPa) 7

Permeability (m2 )

Porosity

Friction angle Dilation angle (◦ ) (◦ )

3.5

1.0 × 10−15

0.12

22

2

3.5 1.375 2.64 7.05 7.05 3.5 1

1.0 × 10−15 1.0 × 10−19 3.43 × 10−16 1.52 × 10−13 1.52 × 10−13 1.0 × 10−15 Host rock permeability Host rock permeability

0.13 0.05 0.12 0.13 0.13 0.12 Host rock porosity Host rock porosity

22 22 11.25 31.87 31.87 22 15

2 2 2 2 2 2 3

30

3

Tensile strength (MPa)

14

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

17

Fig. 5. Contours of final vertical displacement (cm) distributions for the model with weak faults and injection rates of 1 (a), 3 (b) and 5 (c) mt/y at the end of the 20 year injection period. Images in the left column are whole-model 3D plots, and in the right column are the top surface view plots. Modelled uplift increases with injection rate and the distribution with respect to the injection well is slightly different at the lowest injection rate. Note that the colour scales are different.

The properties of faults are chosen based on the aim of simulating weak fault and strong fault scenarios (Table 1). The weak fault scenario with a low friction angle (15◦ ), cohesion of 2 MPa and tensile strength of 1 MPa gives strengths and stiffnesses lower than all the host rocks. In contrast, the strength and stiffness properties for the strong fault are higher than the stratigraphic units of the model (e.g. sandstone units) although friction is similar. While these modelled fault strengths are generally higher than those measured for fault rocks (5–20 MPa, Dewhurst and Jones, 2002; Dewhurst et al., 2002), the point of such high strengths was to investigate whether intact rock in the model outside the fault zones would fail in shear or dilate. In addition, there are reports that fault strength (i.e. shear stress withstood in fault zones) can exceed 100 MPa in strong fault scenarios (e.g. Scholz, 2000) and hydrothermally cemented and healed faults (e.g. Tenthorey et al., 2003). The porosity and permeability of the faults are simulated as being identical to the host rocks, but the code does allow the increase of permeability in the case of fault reactivation and strain localisation. However, it is necessary to note here in advance that the increase of permeability does not occur in these models since faults do not fail or reactivate under the simulated injection conditions as will be shown below. The faults were allowed to transmit

fluids across them with similar permeabilities to the host rock as structural studies suggest the large faults in situ may allow crossfault fluid transmission at reservoir level (Langhi et al., 2013). It was therefore, assumed for ease of modelling that all faults can transmit fluid across them where the reservoir is juxtaposed against reservoir. The present models only simulate a single fluid phase of water for simplicity. As such, the injection of CO2 is actually simulated by the injection of water and we analyse the resultant pressure distributions and flow gradients. The fluid density of 1000 kg/m3 and fluid bulk modulus of 1 × 109 Pa (i.e. approximately intermediate between water and CO2 bulk moduli) are assumed for the fluids in the model. The model essentially simulates one phase flow only. While two phase flow would be expected to reduce permeability in the CO2 plume, the pressure field generated would result in the far field undergoing one phase flow as pore fluid is pushed in front of the advancing CO2 plume. This is a simplifying assumption to ease the modelling process, although potentially the pore pressure distribution might change if the relative permeability issue is significant. However, the experimental study of Levine et al. (2014) showed that the relative permeability of CO2 is clustered around 0.4, suggesting the effective permeability of a water–CO2

18

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

Fig. 6. Comparison of top-surface vertical displacement contours (a) and surface uplift development versus time curves for three locations (b) between the weak fault model and the strong fault model with an injection rate of 5 mt/y. The strong fault models shows slightly less uplift and most modelled uplift is in the first 5 years. Note that the colour scales are different.

system would only be reduced approximately by half from the actual porous rock permeability, less than one order of magnitude. Therefore, the potential pore pressure response to reduced effective permeability would be minimal. In addition, the fluid simulated in the model is more compressible than water and as such this would reduce potential increases in pore fluid pressure coming from other sources. 3.3. Model boundary condition and simulation scenarios “Container”-like static mechanical boundary conditions are adopted in the model. The top of the model is simulated as a free surface, while the four vertical edges of the models are not allowed to move in the direction normal to the edges but are free to move in any direction within the edge planes. The base of the model is not allowed to move in the vertical direction but is free to move in other directions. Gravity is applied in all simulations. Hydrostatic pore fluid pressure gradients are initialised for the model based on in situ conditions. Pore pressures on all the model edges are fixed to simulate permeable or open fluid flow boundary conditions (i.e. allow the exchange of fluids between the internal areas of the model with remote areas outside of the model in the basin). Pore pressures in the internal regions of the model are allowed to change freely, in response to the fluid injection. CO2 injection is simulated by injecting fluid at a constant rate approximately at the Harvey-1 well location at depths of about 2340–2440 m (Fig. 4a). Five injection rate scenarios are investigated, including 1–5 million tonnes per year (mt/y), respectively, which sufficiently cover the range of expected injection rates (Stalker et al., 2013). More specifically, these rates are represented by the injection of 1–5 million cubic metres of fluids per year, respectively, with the specification of fluid density = 1000 kg/m3 . All the simulations are run for an injection period of 20 years.

4. Model results In this study, significant attention is placed on determining the potential excessive mechanical deformation as reflected by ground surface uplift as the result of CO2 injection. Fig. 5 shows the patterns of vertical displacement at the end of the simulated injection period of 20 years for three of the five models with weak faults, with the CO2 injection rates of 1, 3 and 5 mt/y, respectively. We note that some small smooth uplifts have occurred as a result of CO2 injection. In the models assuming weak faults, the maximum surface uplifts (i.e. the highest single-point uplift) are ∼0.46, ∼0.89, ∼1.32, ∼1.75 and ∼2.18 cm for the injection rates of 1–5 mt/y, respectively, at the end of the injection period. The corresponding average surface uplifts over an area of approximately 2.5 km radius around the injection site are ∼0.42, ∼0.76, ∼1.12, ∼1.47 and ∼1.84 cm, respectively, for the five injection rates. It is also noted that for the model with the lowest injection rate of 1 mt/y (Fig. 5a), the maximum uplift location at the surface does not coincide vertically with the location of injection at depth that is between the F5 and F9 faults (see Fig. 4a), but appears to be associated with the structural location of the fault graben bounded by the F5 and F10 faults. With the increase of injection rate to 5 mt/y (Fig. 5b and c), the maximum uplift location at the surface is shifted westwards to approximately coincide vertically with the location of injection at depth. Surface uplift decreases clearly and gradually with the increase of distance away from the injection site. Fig. 6 compares the patterns of surface uplift between the models with weak and strong faults, respectively, using an injection rate of 5 mt/y. There is only a very minor reduction of the maximum surface uplift from 2.18 cm to 2.02 cm (Fig. 6a) when fault strength parameters are changed from weak to strong (Table 1). In addition, the area showing the maximum surface uplift is much smaller in the strong fault model than in the weak fault model. The evolution

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

19

close to the background hydrostatic pore pressure gradients. The relatively rapid pore pressure decrease in the far field is mainly due to the permeable fluid flow boundary conditions at the model edge, which allow some outflow of water in the reservoir horizon. Examination of the mechanical failure states at the end of the modelled injection period show that stratigraphic units and faults have not experienced any geomechanical failure as the result of injection. Fig. 8a presents the geomechanical failure state of the model with weak faults and the highest injection rate, illustrating that all the units and faults are in an intact elastic state with no plastic shear or tensile failure. This is the model with the highest likelihood of failure. Fig. 8b shows the Mohr circles of 3D stresses for five locations defined in the cross section view of Fig. 8a (labelled 1–5). The Mohr circles for the weak faults, F5 and F9 , are clearly not in touch with the failure envelope (location 1 and 2), and those for strong F5 and F9 , are far away from the failure envelope. This indicates that the faults in both weak and strong fault scenarios do not reactivate under our model conditions and will not in model scenarios with lower injection rates either. The 3D stress Mohr circles for location 3 in the Wonnerup sandstone between the F5 and F9 faults are also far away from the failure envelope, confirming the geomechanical stability of the reservoir unit in the model. Similarly, the Mohr circles for location 4 (Yalgorup member) and 5 (Basal Eneabba Shale) are situated safely away from the failure envelopes, showing that the CO2 containment simulated will not be breached under the model conditions. Fig. 9 presents the instantaneous fluid flow velocities at the end of the injection period on an E–W cross section along a central transect of the model (see Figs.1b and 8a). The highest flow velocity range of 10–13.1 m/y is localised immediately around the injection site (Fig. 9a), with velocities rapidly decreasing away from the site. By truncating these high flow velocity values (Fig. 9b), it is possible to show that the low to intermediate flow velocities of 0.6–1 m/y are entirely confined to the compartment of the Wonnerup sandstone between F5 and F9 faults, with flow velocity vectors diverging from the injection site. Flow velocities in far-field areas and in the stratigraphic horizons above the Wonnerup sandstone are small or negligible. There is no visible upward flow discharge along the faults as the faults do not reactivate under the tested scenarios.

5. Discussion Fig. 7. Pore fluid pressure patterns for the −2360 m depth level (injection horizon) at the end of the 20 year injection period for three injection rate scenarios. (a) 1 mt/y, (b) 3 mt/y and (c) 5 mt/y. Background hydrostatic pressure is dark blue. Flow is generally radial. In the high injection rate model, the pressures bank up against the fault due to a narrowing of the flow path across the fault (see Fig. 9). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

of surface uplift with time for three locations around the projected area of injection on the surface (Fig. 6b) shows similar patterns between the weak and strong fault models, with only a small difference in the maximum values. Furthermore, the development of surface uplift with time is heterogeneous and nonlinear. The rate of increase in uplift is significantly greater in the early injection phase (first 5 years) than in the later injection phase. The models also show that injection leads to an increase in pore pressures around the injection site. Fig. 7 illustrates pore pressure distributions at the depth horizon of 2360 m (approximate injection level) in the model at the end of the injection period. Regional pore pressure increases due to injection in the models are small. The maximum increase in the model with weak faults and injection rate of 5 mt/y is ∼1 MPa (Fig. 7c). The highest pore pressure increase anomaly is concentrated to a small area immediately adjacent to the injection site (or well). Pore pressure quickly decreases to levels

Potentially-excessive ground surface uplift presents a major concern for any prospective CO2 injection and storage site. For the SW Hub site, the results of the present models with weak faults suggest that the maximum ground surface uplift is less than 0.9 cm under the scenario of injecting 2 mt/y over a period of 20 years. With the increase of injection rates to 5 mt/y, the maximum uplift reaches ∼2.18 cm, tightly localised in a small area around the injection well (Fig. 6a). These uplift magnitudes are relatively small in comparison with the modelled or measured uplift magnitudes from a number of known case studies. Arsyad et al. (2013) carried out a numerical modelling study of potential ground uplift due to CO2 injection into the Ainoura and Berea sandstone formations in the Nagaoka CO2 injection project (Japan), using the TOUGH2–FLAC3D coupling modelling method (Pruess and Garcia, 2002). They reported the maximum ground uplifts of ∼48 cm for the Berea sandstone and ∼18 cm for the Ainoura sandstone under the scenario with an injection rate of ∼1.103 mt/y (35 kg/s) over the period of ∼20 years. Röhmann et al. (2013) presented the results of a 3D hydro-mechanical model for a rather complex geological structure at the prospective Ketzin CO2 injection project (Germany), using FLAC3D. The maximum ground uplift from their analysis is 2.1 cm for the injection scenario based on pore pressure perturbations simulating CO2 injection over the period of 20 years. The

20

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

Fig. 8. (a) Full 3D model and cross section plots of failure state for the model with weak faults and an injection rate of 5 mt/y, i.e. the model most likely to fail. Number 1–5 in the cross section plot give the locations for plots of 3D stresses in the Mohr circle diagrams. White dot gives the location of the injection point. (b) Presentations of 3D effective stresses in the Mohr circle diagrams for weak fault (F9 and F5 ), strong faults (F9 and F5; dashed lines give the failure envelopes if assuming fault cohesion = 5 MPa), Wonnerup and Yalgorup members, and Basal Eneabba Shale, for respective locations 1–5 as labelled. No failure or reactivation is noted and all stress conditions lie below failure envelopes.

only example (to our knowledge) of comparing modelling results to measured ground surface uplift due to CO2 injection is the In Salah site (Algeria) with the injection rates of 0.5–1.0 mt/y. Satellite measurement data (INSAR) revealed a ground uplift of ∼1.5 cm over the first 3 years and the maximum uplift up to ∼2.5 cm (Vasco et al., 2008; Rutqvist et al., 2010). The numerical modelling study for the In Salah site using the TOUGH–FLAC3D coupling approach (Rutqvist et al., 2010) predicted ground uplifts of 1.2–2.0 cm for a period of 3 years. The maximum ground uplift of ∼2.18 cm predicted by the current model with weak faults and the injection rate of 5 mt/y may

reflect a worst-case scenario for the SW Hub site, based on the following considerations: (1) The model simulates a scenario of injecting CO2 into a single well at a rate of 5 mt/y, whereas in reality, injection is more likely to be carried out over multiple wells with a smaller rate (e.g. 1 mt/y per well over 5 wells). Therefore, potential ground uplift induced by these wells with smaller rates and spreading out over an area could be smaller than the case of the single-well injection scenario with the same amount of CO2 injected.

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

21

Fig. 9. Plots of instantaneous fluid flow velocities (m/y) for an E–W and central cross section (A and B in Fig. 1b, also see Fig. 4) in the model with weak faults and an injection rate of 5 mt/y at the end of the 20 year injection period. (a) Plot of the full velocity range. (b) Plot of truncated flow velocities (i.e. only the velocities below 1.0 m/y), so note the difference in scale. Low flow rates are noted and there is a pinch point at fault F10 , where the flow path area narrows and there is no change in fault permeability.

(2) In the present single-phase fluid model, the injection of CO2 is simulated by water. In reality CO2 has smaller compressibility than water, and hence, the corresponding uplift for the same injection rate should be smaller. On the other hand, CO2 due to its lower density than water and hence greater buoyancy may have some adverse impact due to the effect of compressibility. However, the impact of the buoyancy of CO2 , which has not been modelled in this study would be minimal because the injected CO2 is anticipated to migrate in a mostly flat pattern (like a pancake) with a small column thickness (with respect to the total thickness of the top seal and overlying rocks) and run along the underside of the top seal (e.g. Ennis-King and Wu, 2005). (3) The properties for the weak faults in the model (Table 1) represent a very weak fault scenario; faults are often cemented at a regional scale and have cohesions greater than 2 MPa and friction angles greater than 15◦ (e.g. Muhuri et al., 2003). (4) In the present models, the model base is placed within the Willespie unit. As such, the underburden below the reservoir unit (Wonnerup; see Fig. 4) may not be thick enough. Incorporating a much thicker underburden, which will require additional 3D architecture and geomechanical data for this site, could facilitate the partition of deformation and expansion in reservoir into a thicker underburden and hence reduce surface uplift (Tenthorey et al., 2013). The other major concern for prospective CO2 injection and storage sites evaluated here was potential fault reactivation and breach of the containment formations and hence upward fluid or CO2 seepage. Rinaldi et al. (2013) presented a 2D TOUGH–FLAC coupled geomechanical–fluid flow model for CO2 injection into a generic structure of a horizontally-layered reservoir-caprock sequence simulated as elasto–plastic Mohr–Coulomb materials. They showed that injection with a rate of 1261 t/y (0.04 kg/s) might lead to seismic events with fault slip in the range of 2.27–7.53 cm and fluid leakage through caprocks along the fault with both high and low permeability. However, this essentially represents the worst-case scenario since the fault in their model was simulated as a cohesionless weak zone and the 2D generic structure of the model may also limit the distribution and dispersion of pressures in space

as remarked by the authors. We also carried out 2D simulation of the same architecture at the SW Hub which showed greater surface uplift than 3D models, suggesting 3D models improve on 2D due to more realistic freedom for out of plane movements and 3D pore pressure dissipation. For the SW Hub site, the present 3D models demonstrate that even the scenario with a 5 mt/y single well injection rate and weak faults do not show any fault reactivation or breach of the containment formations, under the given geomechanical and hydrological properties based on in situ data, laboratory measurements and the literature. This can be clearly seen in the model 3D effective stress Mohr-circle patterns (Fig. 8). Effective stresses in the strong faults, the reservoir rock (Wonnerup) and containment horizons (Yalgorup and Basal Eneabba Shale) sit at a safe distance away from the failure envelopes. The stresses for weak faults sit close to (but do not touch) the failure envelopes after injection. In Fig. 8b, the potential failure envelopes for the cases of less strong faults (i.e. fault cohesion = 5 MPa) are inferred and plotted, and the effective 3D stresses in the faults still sit well below failure envelopes. Model fluid flow velocity distributions (Fig. 9) show that flow velocities in the far-field area near model edges are small but non-zero. This suggests that CO2 injection into the reservoir can result in some displacement of water along the reservoir horizon laterally away from the injection site in the basin and hence reduce the likelihood and intensity of significant overpressure development. It should be stressed that the current model results are preliminary with the following limitations. There are still high uncertainties in sub-surface geological structures, particularly for the area east of the Harvey-1 well (Figs. 1 and 4). New geophysical work, in the form of a large-scale 3D seismic survey, has recently been completed for the area and should lead to significant refinement of the structural framework developed by Langhi et al. (2013). There are also uncertainties in the rock and fault geomechanical properties used in the models, although these are tertiary effects in comparison to model geometry and boundary conditions. In addition, cross-fault transmissibility scenarios were not investigated (i.e. fault permeability was the same as host rock permeability) as there were too many uncertainties in estimation of shale gouge ratio (Yielding et al., 1997) from current subsurface data. In this case, faults were allowed to transmit fluids, based on the results

22

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23

from the study of Langhi et al. (2013). In future work where more subsurface data is available and potential cores through faults taken (Sharma, personal communication, 2014), variable fault transmissibility scenarios can be investigated with more rigour. Some of these uncertainties will be addressed by sample acquisition during the next drilling phase anticipated to commence in 2015. Future laboratory experimental work on the geomechanical and hydrogeological parameters for the rocks in the SW Hub will reduce some of the uncertainties (Delle Piane et al., 2013), although they will still be subject to scale issues. The improved availability of data in future may also make it possible to analyse detailed structures in fault zones, in combination with some dynamic fault analyses, similar to the work of Tenthorey et al. (2014) for the Otway Basin CO2 storage site. Finally, the present 3D models are constructed such that stratigraphic and fault structures are perpendicular to the modelled E–W cross section (Fig. 4) in the geological 3D framework (Fig. 3) and as such do not allow stress variations with changes in fault strike. Future modelling efforts will concentrate on the construction of full 3D models from planned 3D seismic data, thus, simulating more realistic regional 3D structures. New information on geological structures and rock properties will be incorporated and a wider range of reservoir and fault permeability change scenarios will be investigated.

6. Conclusions Key conclusions from this modelling study are:

(1) Modelled CO2 injection at the SW Hub results in relatively small ground uplift. In the models with weak faults, average surface uplifts over an area of approximately 5 km diameter around the injection site are ∼0.42, ∼0.76, ∼1.12, ∼1.47 and ∼1.84 cm for the injection rates of 1–5 mt/y at the end of the 20 year injection period, respectively. The corresponding maximum uplifts are ∼0.46, ∼0.89, ∼1.32, ∼1.75 and ∼2.18 cm, respectively. Uplifts are marginally smaller when assuming strong faults. (2) The simulated CO2 injection with a single well of up to 5 mt/y does not lead to the reactivation of faults or breach the containment formations. For this injection rate, Mohr circle plots show that weak faults have stresses close to the failure envelope, but those for other stratigraphic units and strong faults are far away from failure envelopes. (3) The high flow velocities caused by injection are localised immediately around the injection site, with rapid velocity decrease away from the site. Flow patterns suggest that ground waters could be displaced laterally in the Wonnerup sandstone horizon, with flow velocity vectors divergent from the injection site. However, flow velocities in far-field areas and in the stratigraphic horizons above the Wonnerup sandstone are very small or negligible. There is no visible upward flow discharge along the faults and limited pore pressure increase.

Acknowledgements The authors would like to thank CSIRO and the National Geosequestration Laboratory (NGL) for supporting this work. NGL, funded by the Australian Government, is a collaboration between CSIRO, the University of Western Australia and Curtin University. Thanks also go to the Western Australian Department of Mines and Petroleum for allowing data from the SW Hub to be shared and used. We also thank the editor and two anonymous reviewers for suggestions that improved the original manuscript.

References Arsyad, A., Mitani, Y., Babadagli, T., 2013. Comparative assessment of potential ground uplift induced by injection of CO2 into Ainoura, and Berea sandstone formations. Procedia Earth Planet. Sci. 6, 278–286. Aruffo, C.M., Rodriguez-Herrera, A., Tenthorey, E., Krzikalla, F., Minton, J., Henk, A., 2014. Geomechanical modelling to assess fault integrity at the CO2CRC Otway project, Australia. Aust. J. Earth Sci. 61 (7), 987–1001, http://dx.doi.org/10.1080/08120099.2014.958876. Bradshaw, J., Bradshaw, B., Spencer, L., Mackie, V., 2000. GEODISC project 1—regional analysis. Stage 2 Basins—Perth Basin, Western Australia. CO2CRC report number RPT06-0162, 41 pp. Bradshaw, J., Rigg, A., 2001. The geodisc program: research into geological sequestration of CO2 in Australia. Environ. Geosci. 8, 166–176. Causebrook R., Dance T., Bale K., 2006. Southern Perth Basin Site investigation and geological model for storage of carbon dioxide. CO2CRC, report RPT06-0162. Chang, C., Zoback, M.D., Khaksar, A., 2006. Empirical relations between rock strength and physical properties in sedimentary rocks. J. Pet. Sci. Eng. 51, 223–237. Crostella, A., Backhouse, J., 2000. Geology and petroleum exploration of the Central and Southern Perth Basin, Western Australia. Western Australia Geological Survey, report 57, 85 pp. Cundall, P.A., Board, M., 1988. A microcomputer program for modelling large-strain plasticity problems. In: Swoboda, G. (Ed.), Proceedings of the Sixth International Conference on Numerical methods in Geomechanics: Numerical Methods in Geomechanics 6, 2101–2108. Dance, T., Tyson, S., 2006. Geological reservoir characterisation and modelling of the offshore Perth Basin. Research Program 1.1—Technologies for Assessing Sites for CO2 Storage. CO2CRC report number: RPT06-0156. Delle Piane, C., Olierook, H.K.H., Timms, N.E., Saeedi, A., Esteban, L., Rezaee, R., Mikhaltsevitch, V., Iglauer, S., Lebedev, M., 2013. Facies-based rock properties distribution along the Harvey 1 stratigraphic well. CSIRO report EP133710, pp. 156 http://www.anlecrd.com.au/announcements/media-release-research -project-drills-down-for-results Dewhurst, D.N., Jones, M.R., 2002. Geomechanical, microstructural and petrophysical evolution in experimentally reactivated cataclasites: application to fault seal prediction. AAPG Bull. 86, 1383–1405. Dewhurst, D.N., Jones, M.R., Hills, R.R., Mildren, S.D., 2002. Microstructural and geomechanical characterisation of fault rocks from the Carnrvon and Otway Basins. Aust. Petrol. Prod. Explor. Assoc. J. 42, 167–186. Dewhurst, D.N., Hennig, A., 2003. Geomechanical properties related to top seal leakage in the Carnarvon Basin, Northwest Shelf, Australia. Pet. Geosci. 9, 255–263. Ennis-King, J., Wu, G., 2005. Simulation of geological storage of carbon dioxide in the onshore Perth Basin. CSIRO Petroleum, CO2CRC report number RPT05-0087. Fang, H.Y., Daniels, J., 2006. Introductory Geotechnical Engineering. Taylor & Francis, New York, 546 pp. Griffiths, C.M., Seyedmehdi, Z., Salles, T., Dyt C., 2012. Stratigraphic forward modelling for South West Collie Hub phase one–static model. ANLEC R&D project 7-0411-0129; CSIRO report number EP13068. Hennig, A., Otto, C., 2005. A hydrodynamic characterisation of the offshore Vlaming Sub-basin. Research Program 1.2 Technologies for Assessing Sites for CO2 Storage. CO2CRC report number: RPT05-0223. Hillis, R.R., Reynolds, S., 2000. The Australian stress map. J. Geol. Soc. London 157, 915–921. Itasca, 2005. FLAC3D Fast Lagrangian Analysis of Continua in 3 Dimensions, User Manual, Version 3.1. Itasca Consulting Group, Inc., Minneapolis. Langhi, L., Ciftci, B., Strand, J., 2013. Fault seal first-order analysis–SW Hub. CSIRO report EP13879, pp. 50 http://www.anlecrd.com.au/announcements/media-release-research-project -drills-down-for-results Langhi, L., Zhang, Y., Gartrell, A., Underschultz, J.R., Dewhurst, D.N., 2010. Evaluating hydrocarbon trap integrity during fault reactivation using geomechanical 3D modelling: an example from the Timor Sea Australia. AAPG Bull. 94, 567–591. Levine, J.S., Goldberg, D.S., Lackner, K.S., Matter, J.M., Supp, M.G., Ramakrishnan, T.S., 2014. Relative permeability experiments of carbon dioxide displacing brine and their implications for carbon sequestration. Environ. Sci. Technol. 48, 811–818. Mandl, G., 1988. Mechanics of Tectonic Faulting. Elsevier, Amsterdam, 407 pp. Mbia, E.N., Frykman, P., Nielsen, C.M., Fabricius, I.L., Pickup, G.E., Sorensen, A.T., 2014. Modeling of the pressure propagation due to CO2 injection and the effect of fault permeability in a case study of the Vedsted structure, Northern Denmark. Int. J. Greenhouse Gas Control 28, 1–10. McLellan, J.G., Oliver, N.H.S., Schaubs, P.M., 2004. Fluid flow in extensional environments; numerical modeling with an application to Hamersley iron ores. J. Struct. Geol. 26, 1157–1171. Mory, A.J., Iasky, R.P., 1996. Stratigraphy and structure of the onshore Northern Perth Basin, Western Australia. Western Australia Geological Survey, report 46, 101 pp. Muhuri, S.K., Dewers Jr., T.E., Reches, Z.E., 2003. Interseismic fault strengthening and earthquake-slip instability: friction or cohesion? Geology 31, 881–884. Olierook, H.K.H., Delle Piane, C., Timms, N.E., Esteban, L., Rezaee, R., Mory Arthur, J., Hancock, L., 2014. Facies-based rock properties characterization for CO2

Y. Zhang et al. / International Journal of Greenhouse Gas Control 37 (2015) 12–23 sequestration: GSWA Harvey 1 well, Western Australia. Mar. Pet. Geol. 50, 83–102. Ord, A., 1991. Deformation of rock: a pressure-sensitive, dilatant material. Pure Appl. Geophys. 137, 337–366. Pevzner, R., Lumley, D., Urosevic, M., Gurevich, B., Bona, A., Alajmi, M.A., Shragge, J., Pervukhina, M., Mueller, T., Shulakova, V., 2013. Advanced Geophysical Data Analysis at Harvey-1: Storage Site Characterization and Stability Assessment. ANLEC R&D, project number 7-1111-0198 http://www.anlecrd.com.au/ literature 146939/Advanced geophysical data analysis at Harvey-1 - Curtin University Plumb, R.A., 1994. Influence of composition and texture on the failure properties of rock. In: Proceedings of the Eurock SPE/ISRM Symposium on Rock Mechanics, 1994, Delft, SPE 28022, pp. 13–20. Price, N.J., Cosgrove, J.W., 1990. Analysis of Geological Structures. Cambridge University Press, Cambridge U.K, pp. 507pp. Pruess, K., Garcia, J., 2002. Multiphase flow dynamics during CO2 injection into saline aquifers. Environ. Geol. 42, 282–295. Rasouli, V., Pervukhina, M., Müller, T.M., Pevzner, R., 2013. In-situ stresses in the Southern Perth Basin at the GSWA Harvey-1 well site. Explor. Geophys. 44, 289–298 http://dx.doi.org/10.1071/EG13046 Rinaldi, A.P., Jeanne, P., Rutqvist, J., Cappa, F., Guglielmi, Y., 2013. Effects of fault-zone architecture on earthquake magnitude and gas leakage related to CO2 injection in a multi-layered sedimentary system. Greenhouse Gas Sci. Technol. 4, 99–120. Rinaldi, A.P., Rutqvist, J., Cappa, F., 2014. Geomechanical effects on CO2 leakage through fault zones during large-scale underground injection. Int. J. Greenhouse Gas Control 20, 117–131. Röhmann, L., Tillner, E., Magri, F., Kühn, M., Kempka, T., 2013. Fault reactivation and ground surface uplift assessment at a prospective German CO2 storage site. Energy Procedia 40, 437–446. Rutqvist, J., Vasco, D.W., Myer, L., 2010. Coupled reservoir–geomechanical analysis of CO2 injection and ground deformations at In Salah, Algeria. Int. J. Greenhouse Gas Control 4, 225–230. Sarda, J.P., Kessler, N., Wicquart, E., Hannaford, K., Deflandre, J.P., 1993. Use of porosity as strength indicator for sand production evaluation. In: Proceedings of the SPE Annual Technical Conference and Exhibition, Houston, Tex. 3–6 Oct. 1993, Society of Petroleum Engineers, Richardson, Tex. Vol. 5, pp. 381–388. Scholz, C.H., 2000. Evidence for a strong San Andreas fault. Geology 28, 163–166. Sorjonen-Ward, P., Zhang, Y., Zhao, C., 2002. Numerical modeling of orogenic processes and gold mineralization in the southeastern part of the Yilgarn craton, Western Australia. Aust. J. Earth Sci. 49, 1011–1039. Stalker, L., Varma, S., Van Gent, D., Haworth, J., Sharma, S., 2013. South West Hub: a carbon capture and storage project. Aust. J. Earth Sci. 60, 45–58. Strayer, L.M., Hudleston, P.J., Lorig, L.J., 2001. A numerical models of deformation and fluid-flow in an evolving thrust wedge. Tectonophysics 335, 121–145. Tenthorey, E., Cox, S.F., Todd, H.F., 2003. Evolution of strength recovery and permeability during fluid–rock reaction in experimentally fault zones. Earth Planet. Sci. Lett. 206, 161–172.

23

Tenthorey, E., Vidal-Gilbert, S., Backé, G., Puspitasari, R., Pallikathekathil, Z.J., Maney, B., Dewhurst, D.N., 2013. Modelling the geomechanics of gas storage: a case study from the Iona gas field, Australia. Int. J. Greenhouse Gas Control 13, 138–148. Tenthorey, E., Dance, T., Cinar, Y., Ennis-King, J., Strand, J., 2014. Fault modelling and geomechanical integrity associated with the CO2CRC Otway 2C injection experiment. Int. J. Greenhouse Gas Control 30, 72–85. Van Ruth, P., 2006. Geomechanics: Vlaming sub-basin, Western Australia. CO2CRC report RPT06-0043 http://www.co2crc.com.au/dls/pubs/regional/perth/Appendix9 RPT06-0043pdf Vanicek, I., 2013. The importance of tensile strength in geotechnical engineering. Acta Geotech. Slov. 1, 5–17. Varma, S., Dance, T., Underschultz, J., Langford, R.P., Dodds, K., 2007. Regional study on potential CO2 geosequestration in the Collie Basin and the Perth Basin of Western Australia. Prepared for the Department of Industry and Resources and the Coal Futures Group, May 2007. CO2CRC report number: RPT07-0529. Varma, S., Underschultz, J., Dance, T., Langford, R.P., Esterle, J., Dodds, K., Van Gent, D., 2009. Regional study on potential CO2 geosequestration in the Collie Basin and the Southern Perth Basin of Western Australia. Mar. Pet. Geol. 26, 1255–1273. Vasco, D.W., Ferretti, A., Novali, F., 2008. Reservoir monitoring and characterization using satellite geodetic data: interferometric synthetic aperture radar observations from the Krechba field, Algeria. Geophysics 74, WA113–WA122. Verdon, J.P., Kendall, J.M., Stork, A.L., Chadwick, R.A., White, D.J., Bissell, R.C., 2013. Comparison of geomechanical deformation induced by megatonne-scale CO2 storage at Sleipner, Weyburn, and In Salah. Proc. Natl. Acad. Sci. U.S.A. (PNAS), E2762–E2771 www.pnas.org/cgi/doi/10.1073/pnas.1302156110. Vidal-Gilbert, S., Tenthorey, E., Dewhurst, D., Ennis-King, J., Van Ruth, P., Hillis, R., 2010. Geomechanical analysis of the Naylor field, Otway Basin, Australia: implications for CO2 injection and storage. Int. J. Greenhouse Gas Control 4, 827–839. Yielding, G., Freeman, B., Needham, D.T., 1997. Quantitative fault seal prediction. Am. Assoc. Pet. Geol. Bull. 81, 897–917. Zhang, Y., Gartrell, A., Underschultz, J.R., Dewhurst, D.N., 2009. Numerical modelling of strain localization and fluid flow during extensional fault reactivation: implications for hydrocarbon preservation. J. Struct. Geol. 31, 315–327. Zhang, Y., Schaubs, P.M., Zhao, C., Ord, A., Hobbs, B.E., Barnicoat, A.C., 2008. Fault-related dilation, permeability enhancement, fluid flow and mineral precipitation patterns: numerical models. In: Wibberley, C.A.J., Kurtz, W., Imber, J., Holdsworth, R.E., Collettini, C. (Eds.), The Internal Structure of Fault Zones: Implications for Mechanical and Fluid-Flow Properties, 299. Geological Society, London, pp. 239–255, Special Publications. Zhang, Y., Underschultz, J.R., Gartrell, A., Dewhurst, D.N., Langhi, L., 2011. Effects of regional fluid pressure gradients on strain localisation and fluid flow during extensional fault reactivation. Mar. Pet. Geol. 28, 1703–1713.