METHANE BUBBLE ASCENT WITHIN FINE-GRAINED COHESIVE

1 day ago - Methane (CH4) is a potent greenhouse gas. Its release from aquatic sediments to the water column and potentially to the atmosphere, is a ...
0 downloads 0 Views 922KB Size
Subscriber access provided by Bethel University

Environmental Modeling

METHANE BUBBLE ASCENT WITHIN FINE-GRAINED COHESIVE AQUATIC SEDIMENTS: DYNAMICS AND CONTROLLING FACTORS Shahrazad Trboush Sirhan, Regina Katsman, and Michael Lazar Environ. Sci. Technol., Just Accepted Manuscript • DOI: 10.1021/acs.est.8b06848 • Publication Date (Web): 01 May 2019 Downloaded from http://pubs.acs.org on May 2, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 30

Environmental Science & Technology

1

METHANE BUBBLE ASCENT WITHIN FINE-GRAINED

2

COHESIVE AQUATIC SEDIMENTS: DYNAMICS AND

3

CONTROLLING FACTORS Shahrazad Trboush Sirhan, Regina 𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊𝐊∗, Michael Lazar

4 5

The Dr. Moses Strauss Department of Marine Geosciences, The University of Haifa, 199

6

Aba Khoushy Ave., Haifa, Mount Carmel, 3498838, Israel.

7 8

* Corresponding author information:

9

Email: [email protected]

10

Tel.: +972 4 8288979.

11

Fax: +972 4 8288267.

12 13 14 15 16 17

Keywords: bubble, methane, ebullition, muddy aquatic sediment, solute transport, fracture mechanics

18 19

1 ACS Paragon Plus Environment

Environmental Science & Technology

Page 2 of 30

1

ABSTRACT

2

Methane (CH4) is a potent greenhouse gas. Its release from aquatic sediments to the water

3

column and potentially to the atmosphere, is a subject of a great concern. A coupled macroscopic

4

single-bubble mechanical/reaction-transport numerical model was used to explore the ascent of a

5

mature CH4 bubble towards the seafloor in muddy aquatic sediment. Two bubble ascent scenarios

6

were demonstrated: stable and dynamic. For small effective overburden loads (≤ 11 kPa), stable

7

ascent is followed by dynamic ascent (which has not been previously demonstrated to the best of

8

the our knowledge). This ultimately leads to the bubble being released to the water column. Higher

9

effective overburden loads induce only stable bubble ascent, which stops at the gas horizon

10

frequently observed below the seafloor. The depth of the gas horizon increases, while bubble rise

11

velocity decreases with an increase in the overburden load. It is shown that the bubble migration

12

scenario is managed predominantly by inner bubble pressure, which defines a bubble solute

13

exchange with ambient porewaters. Predicting a bubble ascent scenario in muddy sediment will

14

further allow estimation of CH4 emission to the atmosphere and evaluation of changes in the

15

effective mechanical properties of aquatic sediment due to the ascending bubbles.

16 17

2 ACS Paragon Plus Environment

Page 3 of 30

Environmental Science & Technology

1

 INTRODUCTION

2

Methane (CH4) bubbles are abundantly found within aquatic sediment around the world1,2.

3

They are a source of major concern due to their contribution to global warming3,4 and to sediment

4

destabilization. Bubble transport is the significant release mechanism of CH4 from freshwater lakes

5

and reservoirs5,6, wetlands7, and peatlands8.

6

Migration of bubbles significantly affects aquatic environments. Intensive bubble streams

7

induce upwelling flows9, which irrigate sediments, facilitating an exchange of solutes with the

8

overlying porewater10. Bubbles transported in cohesive muddy sediment create sub-vertical

9

fractures11,12 providing lower-resistance conduits for the migration of other bubbles13-17. This can

10

destabilize sediment structures18 and even result in slope failure19,20. In addition, it is suggested that

11

in contrast to dissolved CH4, rising bubbles can bypass processes of oxidation in the upper sediment

12

layers due to their rapid ascent21.

13

Bubble growth in fine-grained cohesive aquatic sediments has been addressed by numerous

14

studies. It has been shown to be managed by a system of coupled mechanical/physical/bio-chemical

15

processes22. Experiments have been performed to characterize bubble growth in muddy

16

sediment13,23, in gelatine (a muddy sediment surrogate)23,24 and in sludge11,25, as well as coupled

17

mechanical and reaction-transport modeling26-30. These studies have described the mechanical

18

interaction of a bubble with surrounding media by Linear Elastic Fracture Mechanics (LEFM), a

19

theory applicable for fine-grained muds consolidated to high yield stress with low fracture

20

toughness31,33. Theory predicts the shape and size of a bubble growing in these media and sheds

21

light on its growth dynamics. Specifically, it has been shown that bubble growth in fine-grained

22

sediments is governed by cavity expansion due to the large capillary-entry pressure precluding gas

23

from entering pore throats11,33,34. In contrast, in coarse-grained sediments, this process occurs by

24

capillary invasion through the sediment framework. Van Kessel and van Kesteren11 and Johnson et

25

al.33 demonstrated experimentally and theoretically that a bubble initially grows quickly to entirely 3 ACS Paragon Plus Environment

Environmental Science & Technology

Page 4 of 30

1

fill the pore space it occupies in fine-grained media and then pushes aside the surrounding sediment

2

matrix. This can initiate a fracture in the stronger sediment (with shear strength on the order of kilo-

3

Pascal), where resistance to expansion is higher.

4

Anderson et al.35 found that bubbles formed in different aquatic sediments exhibit different

5

shapes and sizes. Large isolated gas bubbles in muds, displace both pore water and solid grains and

6

are often non-spherical. Best et al.36 demonstrated that in contrast to small (2 mm diameter) sub-

7

spherical bubbles in silty sands, most bubbles observed in clayey silts appeared as high aspect ratio

8

cavities, up to ~40 mm long, with their longest axes aligned in the sub-vertical plane. Gas injected

9

into a sample of muddy sediment formed flat or cornflake-shaped bubbles14, while a bubble in

10

gelatine resembled an elongated "inverted teardrop"24 with a height to width ratio of ~1 and a much

11

smaller bubble thickness, thus resembling tabular or penny-shaped tensional (Mode I) fractures.

12

In contrast to the bubble growth stage, to date mechanical and reaction-transport controls on

13

bubble ascent towards the sediment-water interface have only been explored separately, despite

14

their coupled nature. For example, diagenetic reaction-transport models37,38 and two-phase models39

15

with no solid mechanics component, allowed for the reproduction of large seasonal shifts in the

16

depth of the sulfate-methane transition zone that explained the temporal position of the gas horizon

17

observed in the field40. Alternatively, Algar et al.41 found that the speed of fracture-driven bubble

18

accent depends on the viscoelastic response of surrounding muddy sediments to stress. However,

19

their model neglected solute exchange between the bubble and ambient pore waters, which may

20

potentially force bubble dissolution42 and even prevent its release into the water column.

21

Despite their importance, bubble rise dynamics and conditions favoring CH4 bubble release

22

into the water column, remain largely unclear. In this study, an approach developed in Katsman et

23

al.29 and Katsman30 is applied, where the fracture-driven growth of a buoyant bubble coupled to

24

solute exchange with surrounding pore-waters, was modeled. For the first time, stable and dynamic

25

mature bubble ascent scenarios are demonstrated and analyzed in the context of possible bubble 4 ACS Paragon Plus Environment

Page 5 of 30

Environmental Science & Technology

1

release to the water column. In contrast to previous studies, the important role of overburden load

2

in bubble solute exchange with surrounding pore waters is specified as well as its effect on the

3

depth of the gas horizon in muddy aquatic sediment.

4 5

 METHODS

6

Model Formulation

7

A single bubble model was developed to simulate buoyancy-driven bubble ascent within

8

strong fine-grained cohesive sediment12 accompanied by fracturing. A sub-vertical bubble (Mode

9

I crack) was modeled, consistent with those observed in some natural sediments and their

10

surrogates13,14,23,24. The model is explained in detail in Katsman et al.29 and Katsman30, where it was

11

used to examine bubble growth from an initial penny-shaped seed to its mature 3D configuration

12

resembling an inverted tear drop with a closed tail in cross-section, just prior to its ascent

13

(Katsman30). The advantage of this model lies in a two-way coupling of a solute transport

14

component37, with a Linear Elastic Fracture Mechanics (LEFM) component27,28. The latter was used

15

to model elastic bubble expansion and discrete differential fracturing of muddy sediments at the

16

evolving bubble front. Table S1 in the Supporting Information presents the equations solved in the

17

model.

18

Discrete differential fracturing modeled over the bubble (Mode I crack) front (Fig.S1 in the

19

Supporting Information) allowed accurate tracking of the evolution of its shape, size, and surface

20

area29,30 required for analysis of bubble rise. It was implemented using principles of LEFM, which

21

postulates that behavior of cracks is defined solely by distribution of the stress intensity factor

22

(SIF), in this case Mode I SIF, 𝐾𝐾𝐼𝐼 , over the evolving crack (bubble) front43,44. Direction of crack

23

propagation and fracturing increment were found using Strain Energy Density Criterion. The

24

increment of differential fracturing at the bubble front is specified in Katsman30. If 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 , which

25

is always attained at the head of an ascending bubble29 (Fig.S1 in the Supporting Information), 5 ACS Paragon Plus Environment

Environmental Science & Technology

1

Page 6 of 30

2

exceeds sediment fracture toughness, 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 (that measures strength of a material at tensile failure),

3

Otherwise it will remain in place.

the bubble will fracture the surrounding sediment along its front and will migrate upwards.

4

Equations described in Table S1 in the Supporting Information were solved at each time step

5

in a fully coupled manner by a direct solver (PARDISO). The numerical model was designed,

6

solved and post-processed in the FEM-based Comsol Multiphysics Simulation Environment, v5.2a.

7 8

Simulation Conditions

9

Simulations were initiated with a predefined, pre-mature bubble with one of its free surfaces

10

located at the symmetry plane of the sediment box (Fig.S1 in the Supporting Information). The

11

bubble was elliptical in its front-view, with an opening (defined by elastic displacement of its

12

surfaces) similar to an inverted tear drop in its cross-section (Katsman30). The tail of the initial

13

bubble was intentionally left open (blue line on the inset of Fig.1a), since its height was a bit smaller

14

than that of a mature bubble (green line on the inset of Fig.1a), as defined in Katsman30. This initial

15

configuration allowed modeling the last stages of bubble growth just prior to the beginning of its

16

rise. The center of the bubble (Z=0 in Fig.S1 in the Supporting Information) was located 17.5 cm

17

below the top surface of the modeled sediment box.

6 ACS Paragon Plus Environment

Page 7 of 30

Environmental Science & Technology

1 2 3 4 5

Fig.1 Initial setting for modeling of an upward migrating bubble under two different ambient dissolved CH4 concentration profiles. Panel (a) shows the initial (t=0 s, blue line) cross-section of the modeled pre-mature bubble (center at Z=0m, head at Z=0.028m, tail at Z=-0.028m) with an open tail and the cross-section of the mature bubble a with closed tail (green line, bubble head at Z≅0.033 m, t=1s for run #1) (enlargement presented in the inset). The difference between the depths of the tails of the pre-mature and mature bubbles is negligible and cannot be recognized in the inset. Panel (b) presents the depth-independent ambient CH4 concentration profile (blue line) and depth-dependent profile (dashed-orange line). In the depth-dependent profile, zero CH4 concentrations were attained above the initial pre-mature bubble head.

6 7 8

The following geochemical, mechanical and numerical input data were used for the numerical simulations: • Sediment bulk density45: 𝜌𝜌 = 1240𝑘𝑘𝑘𝑘 ∙ 𝑚𝑚−3 .

9 10 1

• Young's modulus28,46: 𝐸𝐸 = 5.5 ∙ 105 𝑃𝑃𝑃𝑃, fracture toughness23, 47: 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 = 70𝑃𝑃𝑃𝑃 ∙

12

𝑚𝑚2 , Poisson’s ratio48: ν=0.45.

13

reverences therein).

11

14 15

• Effective porosity of the muddy sediment49: ∅ ≅ 0.2 (see also Sevee50 and • Tortuosity-corrected CH4 diffusion coefficient, D, based on the CH4 molecular

diffusion coefficient51 𝐷𝐷𝑚𝑚 = 10−9 𝑚𝑚2 ∙ 𝑠𝑠 −1. 7 ACS Paragon Plus Environment

Environmental Science & Technology

Page 8 of 30

1

• Two different ambient concentration profiles of dissolved CH4 were modeled

2

throughout the sediment column (Fig.1b): (1) a depth-independent CH4 concentration

3

profile with 𝐶𝐶𝐶𝐶𝐻𝐻4 (𝑎𝑎𝑎𝑎) = 6.5𝑚𝑚𝑚𝑚 = 0.1

4

𝑘𝑘𝑘𝑘 𝑚𝑚3

(at 1m depth, Martens et al.52); and (2) a depth-

dependent (decreasing upward) CH4 concentration profile with maximal 𝐶𝐶𝐶𝐶𝐻𝐻4 (𝑎𝑎𝑎𝑎) =

5

0.1

6

𝑘𝑘𝑘𝑘 𝑚𝑚3

(Fig.1b). The depth-independent profile was associated with a large separation

distance between CH4 production and consumption zones and a small gradient in 𝐶𝐶𝐶𝐶𝐻𝐻4 (𝑎𝑎𝑎𝑎)

7

reduction towards the seafloor, characteristic, for example, of marine sediments42. Source

8

strength, 𝑆𝑆𝐶𝐶𝐶𝐶4 (𝑎𝑎𝑎𝑎) (Eq.S3 in the Supporting Information) was zero in this modeling

9

scenario. The depth-dependent profile was associated with a shorter distance between CH4

10

production and consumption zones and a high gradient in 𝐶𝐶𝐶𝐶𝐻𝐻4 (𝑎𝑎𝑎𝑎) reduction towards the

11

seafloor. This is characteristic of lake sediment38. This profile (Fig.1b) was maintained in

12 13

the sediment box away from the bubble by adjusting the spatial source strength, 𝑆𝑆𝐶𝐶𝐶𝐶4 (𝑎𝑎𝑎𝑎),

14

profile above the initial pre-mature bubble head (Fig.1b).

in the regions of the gradient in CH4 concentrations. 𝐶𝐶𝐶𝐶𝐻𝐻4 (𝑎𝑎𝑎𝑎) = 0 was attained in this

15

• Effective overburden load (combined sediment-water column load) on the top

16 17

surface of the sediment box (Fig.S1 in the Supporting Information) ranged from σ0z =

18

mechanical point of view, soft clay-bearing (cohesive) water-saturated sediment (mud) can

19

be treated as a single (combined) medium22, 53, 54 (see Text S1 in the Supporting Information

20

for more detail), as opposed to separate phases (grains-pore water), as proposed by Terzaghi

21

55

22

particles, 𝐶𝐶𝑠𝑠 , and that of the bulk sediment56, 57, 𝐶𝐶:

23

1.5 𝑘𝑘𝑘𝑘𝑘𝑘 to σ0z ≅ 200 𝑘𝑘𝑘𝑘𝑘𝑘 (see Table 1), based on the following considerations. From a

. Effective stress, 𝜎𝜎 ′ , in muds can be evaluated using the compressibility of solid sediment 𝐶𝐶

𝜎𝜎 ′ = 𝜎𝜎 − �1 − 𝐶𝐶𝑠𝑠 � 𝑢𝑢

(1)

8 ACS Paragon Plus Environment

Page 9 of 30

Environmental Science & Technology

1

where 𝜎𝜎 is total stress and u is pore pressure. In muddy sediments25,34,58, 𝐶𝐶𝑠𝑠 ≅ 𝐶𝐶 ,

which eliminates the effect that pore-water pressure, 𝑢𝑢, has on the effective stress, 𝜎𝜎 ′ :

2

𝜎𝜎 ′ ≈ 𝜎𝜎

3

(2)

As a result, the effective overburden load, 𝜎𝜎𝑧𝑧0 =𝜎𝜎 ′ , on the top surface of the sediment

4 5

box (Fig. S1 in the Supporting Information) combined the weight of both the overlying

6

sediment and water column. Constant overburden load was used in each run (Table 1).

7

• The time step of the calculations was set to equal the time scale of the elastic

8

response in muddy sediment, approximated as ∆𝑡𝑡 = 𝜏𝜏 =

9

𝐻𝐻 𝑈𝑈𝑅𝑅

≅ 0.01𝑠𝑠 (see Katsman et al.29,

for more detail), where H is bubble height and 𝑈𝑈𝑅𝑅 is Raleigh wave velocity in muddy

10

sediment, which represents an upper boundary for dynamic crack propagation under Mode

11

I fracturing44.

12 13

Table 1: Schedule of runs and bubble ascent scenarios observed from the modeling Run

Dissolved

CH4

#

concentration profile

Effective

Interpretation of effective

Observed bubble

overburden

overburden load, 𝜎𝜎𝑧𝑧0 *

propagating

load, 𝜎𝜎𝑧𝑧0 , kPa

scenario

1

Depth-independent

1.5

10 cm s.+3 cm w.

Stable + dynamic

2

Depth-dependent

1.5

10 cm s.+3 cm w.

Stable + dynamic

3

Depth-dependent

11

10 cm s. +1 m w.

Stable + dynamic

4

Depth-dependent

70

10 cm s. +7 m w.

Stable

5

Depth-dependent

99

10 cm s.+10 m w.

Stable

6

Depth-dependent

197

10 cm s.+20 m w.

Stable 9

ACS Paragon Plus Environment

Environmental Science & Technology

1

Page 10 of 30

*s. stands for saturated sediment, w. for water column depth

2 3

 RESULTS

4

The dynamics of ascent of a mature buoyant bubble within fine-grained cohesive aquatic

5

sediments and its controlling factors were analyzed. First, a pre-mature buoyant bubble with a

6

slightly open tail (at the latest stages of its growth, Fig.1a) grew via differential fracturing at its

7

front to attain its mature height, which was pre-defined by the fracture toughness of the muddy

8

sediment (Katsman30). The bubble head migrated upwards from Z=0.028m to Z=0.033m (Fig.1a),

9

whereas its tail migrated downwards by a negligible distance predefined by its small opening

10

(Eq.S2 in Table S1 in the Supporting Information, see also Katsman30). Eventually, compressive

11

stresses at the mature bubble tail forced its closure30, 41, which was preserved until the end of bubble

12

ascent.

13 14

Bubble Propagation Patterns

15

Bubble ascent is controlled by evolution of the stress intensity factor, 𝐾𝐾𝐼𝐼 , at its front (see the

16

Methods Section). Evolution of 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 for a pre-mature buoyant bubble (from the latest stages of

17

its growth, Fig.1a) is shown in Figures 2a and 2b for the depth-independent ambient dissolved CH4

18

concentration scenario and for the depth-dependent one (Fig.1b) respectively, for a small

19

overburden load, 𝜎𝜎𝑧𝑧0 = 1.5 𝑘𝑘𝑘𝑘𝑘𝑘 (runs # 1, 2, Table 1). In both cases, two sequential bubble

20

propagation patterns were observed (Table 1): (1) stable propagation characterized by a saw-tooth

21 22

𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 pattern (at 𝑡𝑡 < 5.1𝑠𝑠 in Fig 2a and 𝑡𝑡 < 12.4 𝑠𝑠 in Fig.2b), followed by (2) dynamic

23

Stable Propagation. During its growth/rise the bubble adsorbs CH4 that induces its

24

propagation characterized by a rising trend line 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 pattern43, 44.

elastic expansion, increases stored mechanical energy, 𝑈𝑈𝑀𝑀 , and subsequently increases 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼

10

ACS Paragon Plus Environment

Page 11 of 30

Environmental Science & Technology

1 2 3 4

(Eq.S2 in the Supporting Information). Once 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 reaches a critical value, 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 ≥ 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 , the bubble experiences an incremental fracturing event. According to LEFM, as long as the

dissipated mechanical energy, 𝑑𝑑𝑑𝑑𝑀𝑀 , is sufficient to generate only the new surface area of the

growing crack increment, dc:

𝑑𝑑𝑈𝑈𝑀𝑀 𝑑𝑑𝑑𝑑

5 6 7 8 9

+

𝑑𝑑𝑈𝑈𝑆𝑆 𝑑𝑑𝑑𝑑

=0

(3)

(where 𝑑𝑑𝑑𝑑𝑆𝑆 is the energy required to create the new surfaces of the growing crack through

breaking of particle-particle bonds43, 59), the fracture event will be followed by a drop in 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼

below 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 . This stage is subsequently followed by another cycle of an increase in 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 due to

persisting adsorption of CH4. Such evolution of 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 is characterized by a saw-tooth ("stable")

10

pattern43,44 (Fig.2). Stable bubble propagation obtained in the modeling (Fig.2), may cease at some

11

distance below the sediment-water interface if solute supply through its surface drops to zero and

12

𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 thus remains below 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 .

13

14 15 16

Fig.2 Mode I Stress Intensity Factor, 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 , at the bubble head as a function of time. Under a small overburden load, 𝜎𝜎𝑧𝑧0 = 1.5𝑘𝑘𝑘𝑘𝑘𝑘 (runs #1 and 2, Table 1), two sequential propagation patterns were observed: (1) Stable propagation and (2) Dynamic propagation. Panel (a) addresses the depth-independent ambient CH4 profile (Fig.1b). Insets (a1) and (a3) show a zoom-in of 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 at the bubble head, while (a2) and (a4) show an enlargement of the corresponding inner bubble pressure, 𝑃𝑃𝑏𝑏 , as a function of time. Panel (b) addresses the depthdependent CH4 ambient profile (Fig.1b), where insets (b1) and (b3) show a zoom of 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 11at the bubble head, and (b2) and (b4) show the corresponding inner pressure of the bubble, 𝑃𝑃𝑏𝑏 . ACS Paragon Plus Environment

Environmental Science & Technology

1 2 3 4

Page 12 of 30

Dynamic Propagation. At a certain stage of bubble ascent (𝑡𝑡 ≅ 5.1𝑠𝑠 in Fig 2a and 𝑡𝑡 ≅

12.4 𝑠𝑠 in Fig.2b) the propagation pattern changes from a stable to a dynamic regime. In contrast to the saw-tooth pattern, dynamic propagation is characterized by a continuous increase in

𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 (when 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 > 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 44), accompanied here by simultaneous solute exchange with the

5

ambient field of dissolved CH4. In this scenario, elastic energy release per increment of crack

6

surface,

7

energy was converted to kinetic energy of crack propagation59,

8 9 10 11 12 13

𝑑𝑑𝑈𝑈𝑀𝑀 𝑑𝑑𝑑𝑑

𝑑𝑑𝑈𝑈𝑀𝑀 , 𝑑𝑑𝑑𝑑

+

exceeded that necessary to form the new surface area of the crack,

𝑑𝑑𝑈𝑈𝑆𝑆 𝑑𝑑𝑑𝑑

+

𝑑𝑑𝑈𝑈𝐾𝐾 𝑑𝑑𝑑𝑑

𝑑𝑑𝑈𝑈𝐾𝐾 𝑑𝑑𝑑𝑑

:

=0

𝑑𝑑𝑈𝑈𝑆𝑆 . 𝑑𝑑𝑑𝑑

Excess

(4)

which accelerated bubble rise. According to LEFM theory, crack propagation may become dynamic if one of the following conditions is met43: (1)

If the crack reaches a critical length. A switch to dynamic propagation may be

realized even under fixed external loading conditions.

14

(2)

15

In the simulations carried out here, the height of the open part of the propagating crack did

16

not change significantly during the period of stable bubble propagation, starting from maturity

17 18

(Fig.1a), which was defined solely by a slightly varying 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 ≅ 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 at its head (Fig.2; see also

19

loading on the surface of the bubble (point 2 above). A bubble approaching the sediment-water

20

interface in the dynamic regime with 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 ≫ 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 at its head will ultimately be released into the

21

If the applied loading is subjected to a rapid time-dependent variation.

Katsman30). Therefore, initiation of dynamic fracturing was probably caused by rapid variations in

water column, as discussed below.

22 23

What Causes a Switch from Stable to Dynamic Propagation? 12 ACS Paragon Plus Environment

Page 13 of 30

Environmental Science & Technology

1

The loading on a bubble’s (crack) surface separating gas from sediment is defined by the

2

difference between the uniform inner bubble pressure, 𝑃𝑃𝑏𝑏 , and local ambient compressive stresses

3 4 5 6 7

within the sediment, 𝜎𝜎𝑦𝑦 , normal to the bubble surface: 𝑃𝑃𝑏𝑏 − 𝜎𝜎𝑦𝑦 (in agreement with a classical

superposition principle of LEFM43). This difference defines an opening at any point on the bubble

surface (Katsman30), and thus 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 at its head (Eq. S2 in Table S1 in the Supporting Information,

see also Katsman30). For instance, larger 𝑃𝑃𝑏𝑏 or smaller 𝜎𝜎𝑦𝑦 generate a bigger bubble opening and

8

thus higher values of 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 . In our case, 𝜎𝜎𝑦𝑦 within the sediment is only a function of the vertical

9

ratio), which in turn increases linearly with depth due to gravity (e.g. Katsman30).

10 11 12 13 14 15 16 17 18 19 20 21 22

local stress component, 𝜎𝜎𝑧𝑧 (e.g. 𝜎𝜎𝑦𝑦 =

𝑣𝑣 1−𝑣𝑣

∙ 𝜎𝜎𝑧𝑧 under uniaxial strain conditions, where 𝑣𝑣 is Poison’s

Alternatively, the change in 𝑃𝑃𝑏𝑏 at fracturing is managed by two processes: (1) a sudden

increase in bubble volume due to differential fracturing at the bubble front that decreases 𝑃𝑃𝑏𝑏 , and (2) solute adsorption that increases 𝑃𝑃𝑏𝑏 . Therefore, the following relation: ∆= ∆𝑃𝑃𝑏𝑏 − ∆𝜎𝜎𝑦𝑦

at

fracturing will define its propagation scenario (here ∆𝑃𝑃𝑏𝑏 is a change in the inner bubble pressure

and ∆𝜎𝜎𝑦𝑦 is the change in the stress component in the sediment at the head of the migrating upward bubble, Fig.S1 in the Supporting Information). If ∆ is negative (𝑃𝑃𝑏𝑏 decreases more than 𝜎𝜎𝑦𝑦 as bubble

ascends, ∆𝑃𝑃𝑏𝑏 < 0, |∆𝑃𝑃𝑏𝑏 | > �∆𝜎𝜎𝑦𝑦 �), then the fracture event is followed by a drop in 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 below

𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 , halting the fracturing phase and starting another phase of elastic expansion. In this case,

bubble propagation is stable. Alternatively, positive ∆ (𝑃𝑃𝑏𝑏 decreases less than 𝜎𝜎𝑦𝑦 as bubble ascends, ∆𝑃𝑃𝑏𝑏 < 0, |∆𝑃𝑃𝑏𝑏 | < �∆𝜎𝜎𝑦𝑦 �) indicates that 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 grows and bubble propagation becomes dynamic

(Fig.2). As a result, because 𝜎𝜎𝑦𝑦 decreases linearly towards the sediment-water interface

(Katsman30), the appearance of dynamic fracturing is only controlled by variations in 𝑃𝑃𝑏𝑏 .

Volume Increase at the Fracturing and its Role in 𝑃𝑃𝑏𝑏 Drop. The drop in 𝑃𝑃𝑏𝑏 is

23

proportional to the increase in bubble volume (Fig.3) at fracturing caused by a displacement of the

24

bubble front by a differential fracturing increment, ∆𝑎𝑎 (Katsman30, Lawn43). At the pre-mature stage 13

ACS Paragon Plus Environment

Environmental Science & Technology

1

Page 14 of 30

2

(Fig.1a), ∆𝑎𝑎 > 0 even at the bubble tail (Fig.S1 in the Supporting Information) and the volume

3

fracturing does not occur). For the pre-mature bubble (Fig.1a), the large increase in volume due to

4

increase at fracturing is larger than that for the mature bubble with ∆𝑎𝑎 = 0 at its tail (where

5

fracturing (Fig.3a) always causes a large inner pressure drop (∆𝑃𝑃𝑏𝑏 < 0, |∆𝑃𝑃𝑏𝑏 | > �∆𝜎𝜎𝑦𝑦 � at the

6

mature bubble with a constant height (Katsman30) propagates mostly by vertical translation, where

7

its widening is relatively small (compared to the pre-mature bubble). The smaller volume increase

8

for the mature bubble (Fig.3b) (combined with high rate of solute adsorption via the bubble surface,

9

discussed below) thus causes a smaller pressure drop (∆𝑃𝑃𝑏𝑏 < 0, |∆𝑃𝑃𝑏𝑏 | < �∆𝜎𝜎𝑦𝑦 � at the bubble head),

10

bubble head), which favors the stable propagation pattern over the dynamic one. Alternatively, a

which may contribute to the appearance of the dynamic propagation pattern.

11 12

Fig.3 Evolution of the volume of the simulated pre-mature bubble (a) and mature bubble (b) (Fig.1a) with time. Fracturing events are marked with red arrows (run #1, Table 1). For the premature bubble, the volume increase at fracturing is larger than that for the mature bubble.

13 14 15

Solute Adsorption Rate and its Role in 𝑃𝑃𝑏𝑏 Drop. In contrast to the role volume increase

16

plays in the drop of Pb at fracturing, higher rates of solute adsorption resulted in a smaller drop in

17

Local solute fluxes over the surface of the bubble are shown in Fig.4a for the depth-

18

Pb , thus advancing the start of dynamic fracturing (see Figs. 2a and 2b).

independent concentration profile, and in Fig.4b for the depth-dependent one (Fig.1b, runs #1 and 14 ACS Paragon Plus Environment

Page 15 of 30

Environmental Science & Technology

1

2, Table 1, respectively), after 9 seconds of bubble propagation in both cases. Colors indicate

2

intensity of the local fluxes over the surface of the bubble.

3

4 5 Fig.4 Plots (a) and (b) present local diffusive fluxes via the surfaces of the bubbles (1/2 of each surface is depicted in front-view) at t=9 s after the beginning of their ascent, under 3 cm water depth and different ambient CH4 concentrations (a: under the depth-independent concertation profile (run #1 in Table 1), b: under the depth-dependent concentration profile) (run #2 in Table 1). Colors indicate intensity of local fluxes. Zero diffusive flux (dark blue) is observed at the closed part of the bubble’s tail. Panels (c) and (d) present the bubble’s opening (the narrow bubble dimension) and propagation distance for the corresponding profiles. Panels (e) and (f) present the corresponding fluxes and bubble opening under a water column of 7 m under the depth-dependent concentration profile (run #4 in Table 1) at t=823 s. See text for explanations. 6

As can be seen in Figure 4a, under depth-independent ambient concentrations the bubble

7

adsorbs solute through the top part of its surface (adjacent to its head), which meets the undisturbed

8

fields of high solute concentration at higher locations during bubble rise. In contrast, for the depth-

9

dependent concentration profile, the bubble adsorbs solute by the bottom part of its surface

10

(adjacent to its tail, Fig.4b), because its head meets lower CH4 concentrations (Fig.1b) or even no 15 ACS Paragon Plus Environment

Environmental Science & Technology

Page 16 of 30

1

dissolved CH4 at its rise. Diffusive flux over the closed part of the bubble tail is zero (shown by a

2

dark blue color in Figs.4a, 4b).

3

In addition, under the depth-independent concentration profile (Fig.4a compared to Fig.4b),

4

the bubble propagated much longer distances relative to its initial position (Fig.1a) as expected.

5

This is mainly due to the earlier appearance of dynamic fracturing, resulting from higher rates of

6

solute adsorption. A larger bubble opening (elastic expansion) and its larger height were also

7 8

observed in this case and correlated with a higher 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝐼𝐼 > 𝐾𝐾𝐼𝐼𝐼𝐼𝐼𝐼𝑖𝑖𝑡𝑡 at the head of the bubble in the

9

The total flux of dissolved CH4 into the bubble via its surface is presented in Fig.5 as a

10

function of time for several time intervals for both ambient concentration profiles (Fig.1b, runs #1

11

and 2, Table 1). In the beginning of bubble ascent, the total diffusive flux grows under both ambient

12

concentration profile cases, whereas under the depth-independent profile it grows faster (Figs. 5a,

13

5b).

dynamic rise scenario (Fig.2, see also Katsman30).

14

15 16 ACS Paragon Plus Environment

Page 17 of 30

Environmental Science & Technology

1 2 3

Fig.5 The total diffusive flux via bubble surface as a function of time depicted for several time intervals. (a) and (b): total fluxes under depth-independent and depth-dependent ambient CH4 concentration profiles (runs #1 and 2, Table 1) for some initial time period (0.2s