RELATING GAS HYDRATE SATURATION TO DEPTH OF SULFATE- METHANE TRANSITION Gaurav Bhatnagar, Walter G. Chapman, George J. Hirasaki Department of Chemical & Biomolecular Engineering Rice University 6100 Main St., Houston, TX, 77005 USA Gerald R. Dickens, Brandon Dugan Department of Earth Science Rice University 6100 Main St., Houston, TX, 77005 USA ABSTRACT Gas hydrate can precipitate in pore space of marine sediment when gas concentrations exceed solubility conditions within a gas hydrate stability zone (GHSZ). Here we present analytical expressions that relate the top of the GHSZ and the amount of gas hydrate within the GHSZ to the depth of the sulfate-methane transition (SMT). The expressions are strictly valid for steady-state systems in which (1) all gas is methane, (2) all methane enters the GHSZ from the base, and (3) no methane escapes the top through seafloor venting. These constraints mean that anaerobic oxidation of methane (AOM) is the only sink of gas, allowing a direct coupling of SMT depth to net methane flux. We also show that a basic gas hydrate saturation profile can be determined from the SMT depth via analytical expressions if site-specific parameters such as sedimentation rate, methane solubility and porosity are known. We evaluate our analytical model at gas hydrate bearing sites along the Cascadia margin where methane is mostly sourced from depth. The analytical expressions provide a fast and convenient method to calculate gas hydrate saturation for a given geologic setting. Keywords: gas hydrate, sulfate-methane transition, modeling, Cascadia Margin Corresponding author: Phone: +1 713 348 5416 Fax +1 713 348 5478 E-mail: gjh@rice.edu NOMENCLATURE j ic Mass fraction of component i in phase j , l m eqbc Methane solubility at the base of GHSZ iD Diffusivity of component i if Normalized flux of component i iF Mass flux of component i g Function denoting integral of porosity term sL SMT depth below seafloor tL Depth of GHSZ below seafloor hL Thickness of gas hydrate layer iM Molecular weight of component i 1 2,Pe Pe Peclet numbers for the two fluxes Q Modified sum of Peclet numbers 1 2,r r Fitting parameters for solubility curve hS Gas hydrate saturation ,f sedU Fluid flux due to sedimentation Proceedings of the 6th International Conference on Gas Hydrates (ICGH 2008), Vancouver, British Columbia, CANADA, July 6-10, 2008. ,f extU Fluid flux due to external source ,f totU Net fluid flux sU Sediment flux z Depth below the seafloor j Density of phase j Porosity Porosity at great depths Subscripts/superscripts/overscript: 0 Value at seafloor g Gas phase h Hydrate phase l , w Water phase or component m Methane component s Sulfate component ~ Denotes normalized or dimensionless value INTRODUCTION Clathrate hydrates of gas, often called gas hydrates, form in pore space of marine sediment along continental margins [1,2]. Their stability depends on temperature, pressure, salinity and gas composition. In general, these conditions restrict gas hydrate occurrence to a finite region below the seafloor, usually referred to as the gas hydrate stability zone (GHSZ) [3]. However, the amount of gas hydrate present within this region (“gas hydrate saturation”) can vary considerably, both globally and locally, because it relates to dynamic inputs and outputs of gas, principally methane, over long (>105 yr) timescales [4,5]. The presence of gas hydrates in marine sediments implies high methane concentrations in pore waters at shallow sub-bottom depths, and a significant methane flux towards the seafloor (Figure 1). This upward methane flux consumes dissolved sulfate, so that most, if not all, seafloor settings with gas hydrate exhibit a relatively shallow and sharp sulfate-methane transition (SMT), a depth interval where pore water SO42- and CH4 concentrations approach zero [6,7,8]. For reasons of mass balance, the depth of the SMT should relate to the uppermost occurrence of gas hydrate (Figure 1) [6]. We have developed a numerical model [9] that has been revised to incorporate a dynamic SMT for systems where methane is supplied from depth [10]. The present analysis differs from previous modeling efforts by developing a complete analytical theory for relating gas hydrate saturation to the depth of the SMT. This approach for quantifying gas hydrate abundance is advantageous because it requires only pore water data from shallow piston cores. Figure 1: (a) Schematic representation of a gas hydrate system showing pore water sulfate and methane concentrations go to zero at some shallow depth below the seafloor. Also shown are different system depths and parameters. (b) Close-up of the sulfate-methane transition (SMT) showing overlap of sulfate and methane profiles. GAS HYDRATE SYSTEMS AND SULFATE DEPLETION In anoxic marine sediments, depletion of pore water sulfate occurs through two main reactions. First, sulfate gets reduced when bacteria utilize solid organic carbon molecules as a substrate [11]. The other important reaction is anaerobic oxidation of methane (AOM), where communities of bacteria and archaea use dissolved methane as follows [12,13,14]: 2 4 4 3 2CH SO HCO HS H O (1) The presence of gas hydrates in shallow sediments implies a significant methane flux towards the seafloor, which can make the second route for sulfate depletion significant [12,15,16,17]. In fact, in regions with even modest upward methane flux, such as Blake Ridge, sulfate reduction can be dominated by AOM. This inference can be supported through careful modeling studies or basic observations [6,15,16,17]. For example, the pore water sulfate gradient can vary by large amounts across regions with methane (e.g., a factor of 16 in the Carolina Rise-Blake Ridge region) despite similar sedimentation rates and total organic carbon (TOC) supply [6,12]. We assume that AOM is the only sulfate sink in our model. This is possible because, as mentioned before, we focus on gas hydrate settings where the TOC content of sediment is low, and where all methane is supplied by deeper sources. MATHEMATICAL MODEL Model framework We first derive a relationship between the depth of the SMT and the upward methane flux using a steady-state mass balance equation for sulfate. This is followed by writing a two-phase methane balance for the system, which links the thickness of the gas hydrate layer and gas hydrate saturation to the methane flux. Finally, by relating the methane flux to sulfate flux at the SMT, we show how the SMT depth is related to thickness of the gas hydrate layer and its saturation. The upward methane flux to shallow sediment depends on the net fluid flux and the methane concentration of rising pore waters. In our modeling, we assume that, at steady-state conditions, gas hydrate extends to the base of GHSZ because of continuous sedimentation. Consequently, pore fluid methane concentration at this depth equals the peak solubility value of the point of three-phase equilibrium [9]. Sulfate mass balance Two assumptions are made in formulating the sulfate mass balance: (1) no sulfate depletion occurs within the sulfate reduction zone (SRZ) due to reduction by solid organic carbon; and, (2) both methane and sulfate react fast enough within the SMT so that their concentrations drop to zero at a single depth. Geochemical data from several gas hydrate settings [6,7,8] indicates that sulfate and methane can co-occur across a horizon several meters thick, suggesting that the species actually react over a finite depth instead of a sharp interface. However, we later normalize all depths by the depth to the base of the GHSZ, which causes the finite SMT transition zone to approach a relatively sharp interface in the dimensionless form. The steady-state sulfate mass balance is: 0 l l s f f s f s cU c D z z , 0 sz L (2) The vertical depth z is set to zero at the seafloor and is positive downwards, following previous work [9]. The mass balance (equation (2)) implies that the mass flux of sulfate, 4SO F , remains constant within the SRZ, and can be rewritten as: 4 l l s f f s f s SO cU c D F z , 0 sz L (3) We recast this relation in dimensionless form. The vertical depth is normalized by tL ( / tz z L ). Consequently, the SMT depth sL is also written in scaled form as /s s tL L L , while sulfate concentration is scaled by 4 0 SOc , its value in standard seawater ( 4 0/l ls s SOc c c ). The net fluid flux ( ,f totU ) can be written as the sum of two components: ,f sedU due to sedimentation- compaction and ,f extU due to upward external flow (Appendix A1). This enables definition of two Peclet numbers that compare each fluid flux to methane diffusion, as follows: , 1 f sed t m U L Pe D , ,2 f ext t m U L Pe D (4) The sulfate balance (equation (3)) can now be rewritten in dimensionless form as (Appendix A1): 4 4 1 2 0 1 1( ) 1 1 l l s s s m SO t f SO m D cPe Pe c D z F L c D ,0 sz L (5) where is the reduced porosity / 1 and is 1 / . The porosity model, assuming hydrostatic pore pressure and equilibrium compaction, and details of non-dimensionalization are given in Appendix A1. To simplify the notation we define the following groups: 1 2 1 ( )Pe Pe Q (6) 4 4 4 0 1 1 SO t SO f SO m F L f c D (7) where Q denotes the modified net fluid flux and 4SO f is a modified sulfate flux. Using these definitions, equation (5) can be written as: 4 1 ll s s s SO cQc D f z (8) where sD = /s mD D . The first boundary condition (B.C.) is applied at the seafloor where the normalized sulfate concentration is equal to unity, while the second is applied at the base of the SRZ (i.e., the SMT), where normalized methane and sulfate concentrations are zero: B.C.: 1lsc at 0z , and (9) B.C.: 0lsc at sz L (10) With these B.C.s, equation (8) can be integrated to give the steady-state sulfate profile [18]: 1 exp ( ) ( ) ( ) 1 exp (0) ( ) s sl s s s Q g z g L D c z Q g g L D , 0 sz L (11) where ( )g z is a function obtained from the integral of the porosity term and is given as: 2 ln( (1 ) (1 ) )( ) 1 zz eg z (12) The following expression for sulfate flux ( 4SO f ) can also be derived as a function of sL and Q using equations (8) and (11): 4 1 exp (0) ( ) SO s s Qf Q g g L D (13) Relationship between sulfate and methane flux We now rewrite the sulfate flux in terms of the methane flux from depth. At the base of the SRZ, the molar fluxes of methane and sulfate are equal, due to the 1:1 stoichiometry of the AOM reaction (equation (1)) [6,12]. Thus, the sulfate mass flux ( 4SO F ) can be written in terms of the methane mass flux ( 4CH F ) from below as follows: 4 4 4 4 SO SO CH CH M F F M , at sz L (14) Substituting equation (14) into equation (7) yields: 4 4 4 4 4 0 1 1 1 SOt SO CH f SO m CH MLf F c D M (15) To simplify the notation, we introduce a dimensionless methane flux 4CH f : 4 4 , 1 1 1 t CH CHl f m eqb m Lf F c D (16) where , l m eqbc is the methane solubility at the base of GHSZ. Using this notation, equation (16) is used to express the dimensionless methane flux in terms of the dimensionless sulfate flux: 4 4 SO CH f f m , where 4 4 4 , 0 l SO m eqb CH SO M c m M c (17) To summarize, using equations (13) and (17), we obtain the following expression between sL and 4CH f : 4 4 / 1 exp (0) ( ) SO CH s s f Q mf m Q g g L D (18) Methane mass balance We now perform similar mass balances on methane and water, and apply them to two distinct spatial domains. The first domain extends from the SMT to the top of gas hydrate, whereas the second domain extends from the top of gas hydrate to the base of the GHSZ (Figure 1). The two-phase (aqueous and hydrate) steady-state methane mass balance, valid from the SMT to the base of the GHSZ, is: 1 0 (1 ) l hs f f m h m h l m h f m UU c S c z cS D z , s tL z L (19) The methane flux invariance can be restated as: 4 (1 ) 1 l l hs m f f m h m h h f m CH U cU c S c S D z F , s tL z L (20) To non-dimensionalize this equation, we utilize the following scalings: , s s f sed UU U , , l l m m l m eqb cc c , , h h m m l m eqb cc c , hh f (21) Using the water mass balance, the methane mass balance (equation (20)) is rewritten in the following form (see Appendix A2 for derivation): 4 1 1 1 1 1 1 l h h ls m h h m w m l m h CH PeUQc S c c c cS f z , 1sL z (22) Dissolved methane zone (SMT to top of gas hydrate) The methane mass balance (equation (22)) is first applied to the region extending from the SMT to the top of the hydrate layer (Figure 1). In normalized form, the thickness of the gas hydrate layer becomes /h h tL L L , while depth to the top of hydrate is (1 hL ). This region ( 1s hL z L ) does not contain any hydrate, so that equation (22) can be simplified by setting 0hS : 4 1 ll m m CH cQc f z , 1s hL z L (23) Methane concentration is zero at the SMT ( sL ), and is equal to the solubility curve at the top of the hydrate layer. Hence, the two boundary conditions for this equation are: B.C.(1): 0lmc at sz L (24) B.C.(2): , 1| h l m m sol Lc c at 1 hz L (25) where , , ,( ) ( ) / l m sol m sol m eqbc z c z c is the normalized methane solubility (mass fraction) in pore water in equilibrium with gas hydrate as a function of the scaled depth ( z ). Analogous to the solution of the sulfate mass balance, equation (23) can be integrated with the above boundary conditions to give the following expressions for methane flux and concentration: 4 , 1| 1 exp (1 ) ( ) h m sol L CH h s Qc f Q g L g L (26) , 1 1 exp ( ) ( ) ( ) | (1 ) 1 exp ( ) h sl m m sol L h s Q g z g L c z c g L Q g L , 1s hL z L (27) Gas hydrate zone (top of hydrate to base of GHSZ) The methane mass balance equation (22) is now applied to the region extending from the top of the gas hydrate layer to the base of the GHSZ (1 1hL z ). In this region, the pore water methane concentration is constrained by the solubility curve, which causes gas hydrate saturation ( hS ) to be the primary dependent variable, which gives the following relation in terms of gas hydrate saturation: 4 1 , , , 1 1( ) 1 1( ) 1 ( ) s m sol h h h h m w m sol h m sol CH PeUQc z S c c c z S c z f , 1 1hL z (28) Several previous simulation results have shown that gas hydrate saturation monotonously increases from zero at the top of the gas hydrate layer to a maximum value at the base of the GHSZ [9,10,18,19,20]. We use this observation to impose the constraint that gas hydrate saturation goes to zero as the top of the hydrate layer is approached through this spatial domain. This condition can be written mathematically as: 0hS as 1 hz L (29) Substituting the above condition in equation (28) gives: 4, , 1( ) ( )m sol m sol CHQc z c z f , as 1 hz L (30) We now have three equations (18, 26, and 30) in terms of four unknowns ( sL , 4CHf , Q and hL ). Hence, by using sL as an input, the other three unknowns can be calculated. Coupled equations for sL and hL In this section, we obtain two non-linear coupled equations in terms of the three variables sL , hL and Q . First, we eliminate 4CH f between equations (18) and (26), which amounts to equating the sulfate flux to the methane flux from depth at the SMT: , 1 |/ 1 exp (1 ) ( )1 exp (0) ( ) hm sol L h s s s QcQ m Q Q g L g Lg g L D (31) Secondly, we equate methane flux in the region containing dissolved methane to the methane flux in the region containing gas hydrate. This helps to eliminate 4CH f between equations (26) and (30), yielding: , 1 , ,1 1 | 1 exp (1 ) ( ) 1| | h h h m sol L h s m sol m solL L Qc Q g L g L Qc c (32) Once sL is known for a particular site, equations (31) and (32) can be solved iteratively (e.g., using a Newton-Raphson or bisection algorithm) to get hL and Q . Gas hydrate saturation profile A major advantage of our formulation is that it gives an analytical expression for the gas hydrate saturation profile (below the top of the hydrate layer) through equation (28). This equation can be rearranged to give the saturation profile as a function of scaled depth z , as follows: 4 , , 1 , , ( ) ( ) 1 1 ( ) 1 ( ) CH m sol m sol h h hs h m w m sol m sol f Qc z c z S z PeU c c c z c z , 1 1hL z (33) RESULTS We first summarize the overall calculation procedure to obtain the results: a. Given ,m solc z and other site-specific parameters ( m , sD , and ), solve coupled equations (31) and (32) to obtain Q and hL . b. Using these values, calculate 4CH f , from any of the three expressions, (18), (26) or (30). c. Substitute into equations (11) and (27) to get the sulfate and methane concentration profiles, respectively. d. By specifying the parameters h , 1Pe , hmc and h wc , equation (33) gives the gas hydrate saturation profile hS z within the GHSZ. Normalized methane solubility curve An important parameter in our formulation is the methane solubility curve within the GHSZ. For sake of demonstration and simplicity, we approximate the solubility curve ,m solc z by an exponential function similar to the form proposed in [21]. We start with this simple two parameter solubility function: 2, 1 r zm solc z re (34) This solubility curve is scaled by methane solubility at the base of the GHSZ, so that its normalized value is equal to unity at 1z . This constraint yields the following relationship between 1r and 2r : 2 1 1 rre 21 rr e (35) which allows us to reduce equation (34) to a single parameter equation: 2 1, r zm solc z e (36) This simple equation (36), with a single fitting parameter, 2r , yields very good fits to solubility curves (Figure 2) obtained through rigorous thermodynamic models (e.g., [9]), for two different seafloor depths, seawater salinity, seafloor temperature of 3°C, and a geothermal gradient of 0.04°C/m. Figure 2: Comparison of normalized methane solubility curves, computed from rigorous thermodynamic models versus those obtained from equation (36). Two different seafloor depths are considered, with the corresponding fitting parameters, 2r , listed in the inset. Effect of sL on the gas hydrate system We now explain the system in terms of the input parameter, sL . The following constant parameter values are assumed for all results shown later: = 6/9, = 9 (which correspond to 0 = 0.7, = 0.1), hmc = 0.134, h = 0.9, 4CHM = 16, 4SOM = 96, seawater sulfate concentration equals 28 mM, and sD = 0.64 [22]. Figure 3: Effect of net fluid flux and variable SMT depths on steady-state sulfate and methane concentration profiles. Specifying sL uniquely constrains the sulfate and methane concentration profiles, as well as the top of the gas hydrate layer. The methane solubility curve corresponds to seafloor depth of 1000m, seafloor temperature of 3°C, and geotherm of 0.04°C/m (Figure 2). Figure 4: Effect of net fluid flux and variable SMT depths on steady-state gas hydrate saturation profiles. Shallow SMT depths indicate higher net methane flux from depth and higher gas hydrate saturation within the GHSZ. 1Pe equals 0.1 for all three cases. Numerical simulation results (crosses) from the model of Bhatnagar et al. [10] match well with the analytical saturation profiles (curves). Figure 3 shows steady-state sulfate and methane concentration profiles, obtained through equations (11) and (27), for three different scaled SMT depths. The solubility curve corresponding to seafloor depth of 1000 mbsl ( 2r = 0.625) is used. Due to co-consumption of sulfate and methane at the SMT, shorter sL indicates higher methane flux from below, thereby leading to a shallower top of the gas hydrate layer (Figure 3). Gas hydrate saturation profiles as a function of scaled SMT depths show higher saturations within the GHSZ with decreasing sL , again due to net increase in methane flux (Figure 4). Increase in the thickness of the hydrate layer with decreasing sL is also evident from the saturation profiles. We further compare steady-state gas hydrate saturation profiles obtained from simulation results (crosses) of Bhatnagar et al. [10], which reveal good agreement between the theory developed in this paper and the numerical formulation. The profiles in Figures 3 and 4 clearly highlight that each distinct value of sL results in a unique profile for dissolved sulfate, methane and gas hydrate saturation. Application to Cascadia Margin sites The Cascadia Margin is an accretionary margin characterized by pervasive upward fluid flow with localized gas venting [7,8]. Results from Ocean Drilling Program (ODP) Leg 204 and Integrated Ocean Drilling Program (IODP) Expedition 311 have given great insight into the complex and heterogeneous gas hydrate distribution at several sites drilled along this margin [7,8,23]. Sites in this region are characterized by relatively high fluid fluxes and low average total organic carbon (TOC) content [8,24], which is indicative of a gas hydrate setting where fluids from depth form the dominant methane source. This makes sites along Cascadia Margin a good location to test our model. We use SMT depths and other data for three Cascadia Margin sites (Table 1) to predict gas hydrate saturations, average saturation and depth to the first occurrence of gas hydrate below the seafloor (Table 2). These sites include ODP Site 889 and IODP Sites U1325 and U1326. Site 889 (ODP Leg 146) has been previously modeled as a gas hydrate system dominated by deeper methane sources [9,20]. Davie and Buffett [20] fit the pore water chloride profile at Site 889 using a coupled numerical model with methane supply from depth. Their results indicate peak hydrate saturation close to 2% at the base of GHSZ and average saturation <1% within the GHSZ [20]. This result agrees favorably with our simulation that shows peak saturation of about 2.7% at the base of GHSZ (Figure 5) and average saturation of 0.6% across the entire GHSZ (Table 2). Hyndman et al. [25] calculated gas hydrate saturation between 25-30% of pore space in the 100 m interval above the base of GHSZ at Site 889 using resistivity log data. However, subsequent calculations using a different set of Archie parameters have revised this estimate to 5-10% in that 100 m interval [26]. Although several parameter uncertainties confront such geochemical and geophysical estimates [8,27], average saturation predicted using our SMT based model concurs with the lower estimates at Site 889. Figure 5: Steady-state gas hydrate saturation profiles computed from scaled SMT depths at Cascadia Margin Sites 889, U1325 and U1326. Scaled SMT depth is highest for Site 889 and lowest for Site U1326, implying higher methane flux and greater gas hydrate saturation at Site U1326 and relatively low methane flux and hydrate saturation at Site 889. For the IODP Expedition 311 sites, drilled along the northern Cascadia Margin, we compare our predictions with average saturations computed from chloride anomalies and resistivity log data (Table 2). Average saturation is calculated from chloride data by assuming a background in situ chloride profile and attributing the relative pore water freshening to gas hydrate dissociation (e.g., [27]). Average saturation is obtained from resistivity data using the Archie equation and parameters given in [8]. Average saturation over the GHSZ at Site U1325 is estimated from resistivity data and chloride anomalies to be 3.7% and 5.3%, respectively. Corresponding estimates from resistivity and chlorinity for Site U1326 are 6.7% and 5.5%, respectively. These values compare favorably with 3.1% and 6.6% average saturation from our SMT based model at Sites U1325 and U1326, respectively (Table 2). In general, we get good first order agreement between average gas hydrate saturations derived using resistivity logs/chloride anomalies and those predicted using our model, although our model consistently predicts lower average saturation at all three sites along Cascadia Margin. For estimates from resistivity logs, a possible explanation for the deviation is that interpretations of resistivity logs depend on knowledge of formation water resistivity and three empirical constants, which are hard to constrain in clay-rich sediments. Moreover, most studies employing transport models (e.g., [9,10,20,28]) predict gas hydrate to first occur well below the seafloor. In contrast, log-based results often predict gas hydrate starting immediately below the seafloor. This will cause saturations from transport models to be lower than those predicted using resistivity log data. Similarly, estimation of hydrate saturation from chloride data is quite sensitive to the choice of baseline curves. Apart from the small deviations between model and chloride/resistivity log predictions, our model gives a good average estimate of gas hydrate saturation. CONCLUSIONS We have developed analytical expressions to estimate gas hydrate saturation from scaled depth of the sulfate-methane transition (SMT) for gas hydrate systems dominated by deep-methane sources. This scaled SMT depth is the ratio of the dimensional depth of the SMT below the seafloor to the depth of the gas hydrate stability zone (GHSZ) below the seafloor. Using simple one- dimensional mass balances for sulfate and methane, we show that net methane flux in such deep-source systems uniquely determines the scaled SMT depth, the thickness of the gas hydrate layer and gas hydrate saturation within the GHSZ. Steady-state results show that as the SMT becomes shallower, methane flux and, consequently, gas hydrate saturation increases. Average saturations over the GHSZ at three Cascadia Margin locations, calculated from our method, are 0.6%, 2.3% and 5.5% for Sites 889, U1325 and U1326, respectively. These values compare favorably with averages computed from resistivity log and chlorinity data for all sites. Hence, our analytical formulation provides a simple and fast technique to constrain gas hydrate saturation in deep-source systems. ACKNOWLEDGEMENTS We acknowledge financial support form the Shell Center for Sustainability, the Kobayashi Graduate Fellowship, and the Department of Energy (DE- FC26-06NT42960). REFERENCES [1] Kvenvolden KA. Gas hydrates: Geological perspective and global change. Reviews of Geophysics 1993; 31:173-187. [2] Sloan ED, Clathrate hydrates of natural gases, Marcel Dekker, New York, 1998. [3] Dickens GR. The potential volume of oceanic methane hydrates with variable external conditions. Organic Geochemistry 2001; 32:1179- 1193. [4] Dickens GR. Rethinking the global carbon cycle with a large, dynamic and microbially mediated gas hydrate capacitor. Earth and Planetary Science Letters 2003; 213:169-183. [5] Buffett BA, Archer D. Global inventory of methane clathrate: Sensitivity to changes in the deep ocean. Earth and Planetary Science Letters 2004; 227:185-199. [6] Borowski WS, Paull CK, Ussler WIII. Global and local variations of interstitial sulfate gradients in deep-water, continental margin sediments: Sensitivity to underlying methane and gas hydrates. Marine Geology 1999; 159:131-154. [7] Tréhu AM, Bohrmann G, Rack FR, Torres ME, et al. (Eds), Proceedings of the Ocean Drilling Program, Initial Reports 2003, vol. 204, Ocean Drilling Program, College Station, TX. [8] Riedel M, Collett TS, Malone MJ, and the Expedition 311 Scientists (Eds). Proceedings of the Integrated Ocean Drilling Program 2006, vol. 311, Integrated Ocean Drilling Program Management International Inc., Washington, DC. [9] Bhatnagar G, Chapman WG, Dickens GR, Dugan B, Hirasaki GJ. Generalization of gas hydrate distribution and saturation in marine sediments by scaling of thermodynamic and transport processes. American Journal of Science 2007; 307:861-900. [10] Bhatnagar G, Chapman WG, Dickens GR, Dugan B, Hirasaki GJ. Sulfate-methane transition as a proxy for average methane hydrate saturation in marine sediments. Geophysical Research Letters 2008; 35:L03611, doi:10.1029/2007GL032500. [11] Berner RA. Early Diagenesis: A Theoretical Approach, Princeton Univ. Press, Princeton, N.J., 1980. [12] Borowski WS, Paull CK, Ussler WIII. Marine pore-water sulfate profiles indicate in situ methane flux from underlying gas hydrate. Geology 1996; 24(7):655-658. [13] Reeburgh WS. Methane consumption in Cariaco trench waters and sediments. Earth and Planetary Science Letters 1976; 28:337-344. [14] Orphan VJ, House CH, Hinrichs K-U, McKeegan KD, DeLong EF. Methane-consuming archaea revealed by directly coupled isotopic and phylogenetic analysis. Science 2001; 293:484-487. [15] Davie MK, Buffett BA. A steady state model for marine hydrate formation: Constraints on methane supply from pore water sulfate profiles. Journal of Geophysical Research 2003; 108:2495, doi:10.1029/2002JB002300. [16] Luff R, Wallman K. Fluid flow, methane fluxes, carbonate precipitation and biogeochemical turnover in gas hydrate-bearing sediments at Hydrate Ridge, Cascadia Margin: Numerical modeling and mass balances. Geochimica et Cosmochimica Acta 2003; 67:3403-3421. [17] Snyder GT, Hiruta A, Matsumoto R, Dickens GR, Tomaru H, Takeuchi R, Komatsubara J, Ishida Y, Yu H. Pore water profiles and authigenic mineralization in shallow marine sediments above the methane-charged system on Umitaka Spur, Japan Sea. Deep-Sea Res. II 2007; 54:1216-1239. [18] Bhatnagar G. Accumulation of gas hydrates in marine sediments, Ph.D. Thesis, Rice University, Houston, TX, February 2008. [19] Davie MK, Buffett BA. A numerical model for the formation of gas hydrate below the seafloor. Journal of Geophysical Research 2001; 106:497-514. [20] Davie MK, Buffett BA. Sources of methane for marine gas hydrate: inferences from a comparison of observations and numerical models. Earth and Planetary Science Letters 2003; 206:51-63. [21] Davie MK, Zatsepina OY, Buffett BA, Methane solubility in marine hydrate environments. Marine Geology 2004; 203:177- 184. [22] Iversen N, Jørgensen BB. Diffusion coefficients of sulfate and methane in marine sediments: Influence of porosity. Geochimica et Cosmochimica Acta 1993; 57:571-578. [23] Tréhu AM, et al. Three-dimensional distribution of gas hydrate beneath southern Hydrate Ridge: Constraints from ODP Leg 204 . Earth and Planetary Science Letters 2004; 222:845-862. [24] Westbrook GK, Carson B, Musgrave RJ, et al. (Eds.), Proceedings of the Ocean Drilling Program, Initial Reports 1994, vol. 146 (Pt. 1), Ocean Drilling Program, College Statiom, TX. [25] Hyndman, RD, Yuan T, Moran K. The concentration of deep sea gas hydrates from downhole electrical resistivity logs and laboratory data. Earth and Planetary Science Letters 1999; 172(1):167-177. [26] Collett, TS. Quantitative well-log analysis of in-situ natural gas hydrates, Ph.D. Thesis, Colorado School of Mines, Golden, CO, 2000. [27] Egeberg PK, Dickens GR. Thermodynamic and pore water halogen constraints on hydrate distribution at ODP Site 997 (Blake Ridge). Chemical Geology 1999; 153:53-79. [28] Xu W, Ruppel C. Predicting the occurrence, distribution, and evolution of methane gas hydrate in porous marine sediments. Journal of Geophysical Research 1999; 104:5081-5096. APPENDIX A1. Non-dimensionalization of Sulfate Mass Balance The net fluid flux in the system ( ,f totU ) results from the combination of fluid flux due to continuous sedimentation and compaction of sediments ( ,f sedU ) and the external fluid flux ( ,f extU ) [9,15,18]: , , ,f tot f sed f extU U U (A1) In terms of Peclet numbers, equation (4), this sum can be written as: , , , 1 2 f tot t f sed t f ext t m m m U L U L U L Pe Pe D D D (A2) Multiplying equation (3) by /t mL D and dividing by 4 0 SOc gives: 4 4 1 2 0 l SOl s s t s m f SO m FD c LPe Pe c D z c D (A3) Porosity loss is modeled by relating it to effective stress and assuming hydrostatic pressure (equilibrium compaction), which yields the following relationship between the reduced porosity and normalized depth [9,18]: 1 ze (A4) where and are reduced porosities defined in terms of the maximum ( 0 ) and minimum ( ) porosities achieved during compaction: 1 , 0 1 (A5) Dividing equation (A3) by (1 ) we obtain the dimensionless form of the sulfate mass balance, which is equation (5) in the main text. A2. Non-dimensionalization of methane mass balance The steady-state water mass balance below the SMT can be written as: 0 1 l hs f f w h w h UU c S c z , s tL z L (A6) Equation (A6) can also be written in terms of the water flux ( 2H O F ) as: 2 , , 1 l hs f f w h w h H O f sed f ext f UU c S c F U U (A7) Due to low methane solubility in water, we assume the mass fraction of water in aqueous phase to be unity. This gives us an expression for the water flux: , , 1 hs h f f sed f ext h w f UU U U S c (A8) Substituting this expression for fluid flux into equation (20), we get: 4 , , 1 1 1 h ls h f sed f ext h w m f l CHhs h m h m h m f f UU U S c c FU cS c S D z (A9) Similar to the sulfate mass balance, we multiply the above equation by /t mL D and divide by , l m eqbc to get the following dimensionless form: 4 1 1 1 2 , 1 1 1 h ls s h w h m l CHh m t h m h h l f m eqb m PeU PeUPe Pe S c c Fc LS c S z c D (A10) Equation (A10) can be divided by (1 ) to express in terms of the reduced porosity and rearranged to get the dimensionless methane balance, equation (22). TABLES Table 1: Site-specific parameters for Cascadia Margin sites * Calculated from thermodynamic model [9] ** Calculated from fitting equation (36) to solubility curves obtained from thermodynamic model [9] a ODP Leg 146 [24] b IODP Expedition 311 [8] c S was not available, hence assumed equal to nearest site U1325. Table 2: Results for Cascadia Margin sites Site 1Pe sL t h L L (m) h GHSZS (calc.) h GHOZS (calc.) h GHSZS (res. log) a h GHSZS (Cl-) b 889 0.068 0.058 113 0.6% 1.1% - <1% U1325 0.11 0.059 62 2.3% 3.1% 3.7% 5.3% U1326 0.11 0.060 36 5.5% 6.6% 6.7% 5.5% a From Archie equation using LWD log data [8] b From fit to chloride data Site S (cm/k.y.) T0 (°C) G (°C/m) D0 (m) sL (m) tL (m) , l m eqbc * m (eq. 17) 2r ** (eq.36) 889a 25 3 0.054 1311 10 225 2.1×10-3 4.4 0.73 U1325b 38.3 3 0.06 2195 5 230 2.5×10-3 5.2 0.93 U1326b 38.3c 3 0.06 1828 2.5 126 2.3×10-3 4.8 0.86