# Investigating the velocity of magmatic intrusions and its relation with rock fracture toughness: insights from laboratory experiments and numerical models

Ayleen Gaete, Francesco Maccaferri, Séverine Furst, Virginie Pinel

## ► To cite this version:

Ayleen Gaete, Francesco Maccaferri, Séverine Furst, Virginie Pinel. Investigating the velocity of magmatic intrusions and its relation with rock fracture toughness: insights from laboratory experiments and numerical models. Geophysical Journal International, 2025, 240 (1), pp.638-651. 10.1093/gji/ggae396. ird-04995410

HAL Id: ird-04995410

https://ird.hal.science/ird-04995410v1

Submitted on 18 Mar 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Investigating the velocity of magmatic intrusions and its relation with rock fracture toughness: insights from laboratory experiments and numerical models

Ayleen Gaete<sup>1,2,3</sup> Francesco Maccaferri<sup>1</sup> Séverine Furst<sup>1,4</sup> and Virginie Pinel<sup>1,5</sup>

Accepted 2024 October 27. Received 2023 October 21; in original form 2024 September 11

## SUMMARY

A key question for those who study magmatic and volcanic processes is: ‘How fast can a magmatic intrusion travel?’ Observations and models indicate ranges between $10^{-2}$ and $1 \text{ m s}^{-1}$ depending on several parameters, including magma buoyancy (or driving pressure), viscosity and rock fracture toughness ($K_c$). However, $K_c$ values are difficult to constrain, as effective values inferred from large magmatic intrusions may be 2–3 orders of magnitude larger than measured values from small laboratory samples. This can be attributed to non-elastic processes that dissipate energy at different rates, depending on factors such as the fracture dimension and fracture propagation velocity. Here, we aim to investigate this aspect and provide a scheme for estimating effective fracture toughness values ($K_{\text{eff}}$) by considering fluid-filled fracture processes across different ranges of propagation velocities. To do so, we combine (i) analogue laboratory experiments involving the propagation of oil- and air-filled cracks within a solidified gelatin block, with (ii) numerical simulations, reproducing the crack shape and velocity and providing an estimate of the energy dissipated by the fluid flow between the crack walls. We show that even at the scale of our experiments, $K_{\text{eff}}$ values exhibit significant variations spanning over an order of magnitude. Over the velocity ranges relative to our two sets of experiments, we identify two empirical relations for an effective, velocity-dependent fracture energy ($\Delta E_f(v)$), showing that when such an empirical relation is implemented into the numerical model, it improves the prediction of velocities and velocity variations. Following a similar procedure and building empirical relations for $\Delta E_f(v)$ or $K_{\text{eff}}(v)$ at the scale of magmatic intrusions would improve predictions on dyke propagation velocities in the crust. In order to do so, a considerable amount of observations on the geometry and propagation velocity of magmatic dykes should be gathered.

**Key words:** Fracture and flow; Numerical modelling; Experimental volcanism; Physics of magma and magma bodies.

## 1. INTRODUCTION

Magma propagation through the Earth’s crust mainly occurs by dyking, which may either feed eruptive vents or stall without reaching the surface. Dyking processes are mainly driven by the interplay between crustal stress orientation, rock mechanics and magma properties. Addressing the dynamic processes involved in the propagation of magmatic dykes has been the focus of several volcanological studies using theoretical, laboratory and observational approaches (e.g. Rivalta et al. 2015).

In particular, determining the propagation velocity of an intrusion, given the characteristic parameters of magma and host rocks, currently represents a problem of fundamental interest for volcanology (e.g. Rivalta et al. 2015; Möri & Lecampion, 2022, 2023; Davis et al. 2023; Furst et al. 2023). This is partially due to the relatively scarce amount of direct observations but also to the theoretical challenges that models have to overcome to fully describe the physics of fluid-filled fracture propagation. From an observational perspective, monitored volcanoes offer the best opportunities. For instance, the hypocentre migration of earthquakes associated with the cracking

<footnote n="1"><sup>1</sup>INGV – Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Vesuviano, via Diocleziano, 328 - 80124, Naples, Italy. E-mail: ayleen.gr@gmail.com (AG); francesco.maccaferri@ingv.it (FM)</footnote>

<footnote n="2"><sup>2</sup>GFZ – German Research Centre for Geosciences, Telegrafenberg 14473, Potsdam, Germany</footnote>

<footnote n="3"><sup>3</sup>OVDAS – Observatorio Volcanológico de Los Andes del Sur, SERNAGEOMIN, Rudecindo Ortega 03850, Temuco, Chile</footnote>

<footnote n="4"><sup>4</sup>GEOMAR Helmholtz – Centre for Ocean Research, Wischhofstraße 1-3, 24148, Kiel, Germany</footnote>

<footnote n="5"><sup>5</sup>Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, Univ. Gustave Eiffel, ISTerre, 38000, Grenoble, France</footnote>
at the tip of an intrusion allows for tracking propagation velocities, which have been found to range from km d⁻¹ to km hr⁻¹ (≈ 10⁻² to 1 m s⁻¹, e.g. Prejean et al. 2003; Sigmundsson et al. 2015, 2024; Lengliné et al. 2021).

Theoretical models based on the lubrication theory describe the velocity of fluid-filled fracture propagation under conditions where viscous forces outweigh fracture resistance. This occurs when the velocity is determined by the fluid viscosity (Spence et al. 1987; Lister 1990; Spence & Turcotte 1990; Rubin 1998; Roper & Lister 2007; Davis et al. 2023; Möri & Lecampion 2023). However, they fail when viscous forces are not the limiting factor for the propagation velocity of an intrusion, and the velocity of fracture growth is instead governed by rock fracture toughness ($K_c$), which represents the resistance of the rocks to fracture (Garagash & Germanovich 2014, 2022; Möri & Lecampion 2022, 2023, and references therein). Magma viscosity and rock fracture toughness values span several orders of magnitude. This suggests that dyking processes may occur at different propagation regimes, ranging from viscous- to fracture-dominated (Rivalta et al. 2015). Recently, new 2-D and 3-D models have been developed to simulate buoyant fluid-filled fracture propagation within these two limits. Möri & Lecampion (2023) addressed the release of a finite-fluid volume within a 3-D planar fracture. The fluid-filled fracture self-propagates when it overcomes the hydrostatic limit, characterized by the buoyancy and the viscosity factors. Importantly, the viscosity factor also depends on the release rate of fluid, which in turn influences the propagation regimes and affects the fracture shape. Noteworthy, they found that even in the case of a finite volume, the propagation regime is likely to switch between the two limits. On the other hand, Furst et al. (2023) implemented a new 2-D numerical approach to estimate the propagation velocity of viscous fluid-filled cracks within a regime between fracture- and viscous-dominated. This model considers both rock fracture toughness and fluid viscosity, resulting in velocities and shapes that depend on both factors.

Complementary to numerical models, analogue laboratory experiments—performed by injections of different fluids into solidified gelatin blocks—are additional powerful tools to investigate magmatic dyke dynamics. These fluid-filled fracture experiments share similar physics with magmatic dykes, as they both propagate as quasi-static fluid-filled cracks, with negligible contribution of inertial terms (Rivalta et al. 2015; Kavanagh et al. 2018). Some of these experiments have shown that the ratio between fracture toughness ($K_c$) and stress intensity factor ($K_I$)—representing the intensity of the stress singularity at the crack tip—correlates with the propagation velocity of fluid-filled fractures (Takada 1990; Heimpel & Olson 1994; Smittarello et al. 2021). These results confirm that non-elastic processes occurring at the propagating tip of the crack determine its propagation velocity, as also noted by Rubin (1993). Additionally, similar analogue experiments have been used in the past to study velocity variations when a fluid-filled crack approaches the surface (Watanabe et al. 2002; Rivalta & Dahm 2006), crosses a layer with different elastic properties (Maccaferri et al. 2010), or enters a heterogeneous stress field (Pinel et al. 2022). Eventually, previous theoretical and empirical results suggest that when the velocity of a fluid-filled fracture is not constrained by fluid viscosity, the fundamental parameter governing the propagation velocity is the fracture toughness of the host rocks (Heimpel & Olson 1994; Garagash & Germanovich 2014; Möri & Lecampion 2023; Furst et al. 2023). Theoretical and experimental studies on mode-I fracture mechanics showed that the process of tensile fracturing—under different loading stress conditions—is governed

by velocity-dependent empirical relations for static and dynamic rock fracture toughness, as well as crack stress intensity factor (Bhat et al. 2012, and references therein). However, these tensile fracturing processes do not share the same driving mechanism with fluid-filled cracks, and importantly, they may occur at very different velocity regimes, with significantly higher velocities resulting in differences between the static and dynamic fracture toughness. Monitoring hydraulic fracture propagation rates in rocks remains challenging (Liu & Lecampion 2022, and references therein), and only a few studies were able to highlight velocity-dependent characteristics of the fracture growth process within the velocity ranges for fluid-filled fracture propagation (Chen et al. 2021, and references therein). Recently, Liu & Lu (2023) investigated the effects of velocity-dependent apparent fracture toughness on hydraulic fracture propagation using a power-law empirical model. They present results for plane strain and radial fractures, including an application to a non-buoyant magmatic intrusion under constant magma supply. They found that both plane strain and radial fractures present decreasing toughness as the fracture grows. Additionally, for the radial fracture, the energy spent in creating new fracture surfaces increases as the fracture extends, eventually causing the propagation to evolve towards a regime dominated by fracture toughness, in contrast to the plane strain case. However, estimates of rock fracture toughness can vary by orders of magnitude depending on the scale of the fracture process. Typically, laboratory estimates on small, centimetre-scale rock samples return fracture toughness values that are 2–3 orders of magnitude smaller than those inferred from large, kilometre-scale, magma-filled dykes exposed in the field (Delaney & Pollard 1981; Olson & Schultz 2011; Zhu et al. 2022). Such different estimates are thought to be related to effective fracture toughness values resulting from non-elastic (velocity- or scale-dependent) processes occurring at the tip of a propagating crack (Delaney et al. 1986; Rubin 1993, 1995; Rivalta et al. 2015).

In this study, we conducted two sets of analogue experiments involving the injection of fluids with extremely different viscosities—silicon oil and air—into several gelatin blocks with very similar elastic properties. We performed multiple injections with different volumes for each fluid, resulting in a wide range of propagation velocities spanning over 4 orders of magnitudes. We used the experiment records to measure the crack geometries and estimate their stress intensity factors. We also used the numerical code developed by Furst et al. (2023) to reproduce the propagation velocity of oil-filled fractures with different volumes and derived their stress intensity factor ($K_I$) accounting for the viscous oil flow between the fracture walls. In fact, during previous analogue experiments of fluid-filled crack propagation, $K_I$ values have often been estimated based on static crack formulae, which neglect the effect of fluid viscosity (e.g. Takada 1990; Heimpel & Olson 1994; Smittarello et al. 2021). We used $K_I$ estimates from oil- and air-filled cracks to compute the gelatin fracture toughness ($K_c$) and eventually showed that these two sets of experiments yield notably different $K_c$ estimates. We interpreted our findings in terms of velocity-dependent effective fracture toughness values ($K_{eff}$), with slower fracture propagation velocities requiring less energy to break a unit surface of gelatin at the crack tip. This is consistent with previous experimental results (Heimpel & Olson 1994) and theoretical considerations (Rubin 1993). Our results suggest that introducing non-elastic, velocity- (or, in general, ‘crack-’) dependent effective fracture toughness is fundamental for theoretical models to predict fluid-filled fracture propagation velocities successfully.
# 2. METHODS

## 2.1. Analogue model framework

We conducted experiments using solidified gelatin as the analogue material for the Earth's crust. Gelatin is transparent and brittle-elastic, and it is broadly used to study magma transport by injecting different fluids (e.g. air, water and silicon oils) to form fluid-filled fractures. These experiments are representative of intrusion mechanisms that involve brittle-elastic properties of the crust, offering a good scaling for different types of magma (e.g. Takada 1990; Heimpel & Olson 1994; Watanabe et al. 2002). To ensure adequate scaling for magmatic dyke propagation studies, gelatin must be prepared at a concentration < 5 wt% and cooled down at a temperature below 10 °C (Kavanagh et al. 2013).

To simulate the mechanics of magmatic dykes at toughness- and viscosity-dominated regimes, we respectively injected air and silicon oil (specifically, KORASINOL M10000 supplied by Obermeier, which has a viscosity 10000 times higher than water) into the gelatin, resulting in the formation of fluid-filled buoyant cracks. Gelatin and silicon oil have similar densities ($\rho \sim 1000 \text{ kg m}^{-3}$ and $\rho_{\text{oil}} \sim 980 \text{ kg m}^{-3}$ at 10 °C, respectively); therefore, we increased the buoyancy of the oil-filled cracks by dissolving salt within the gelatin during preparation. While salt increases the gelatin density, improving the model scaling for buoyant magma (Brizzi et al. 2016), it also affects the gelatin structure, decreasing its rigidity (Smittarello et al. 2021). Furthermore, previous studies showed that—for a given shear modulus—salty gelatin blocks display higher values of fracture toughness with respect to non-salty gelatin (Smittarello et al. 2021).

### 2.1.1. Gelatin preparation and laboratory setup

The experiments were performed in the Volcano Lab at the German Research Centre for Geoscience (GFZ) in Potsdam, Germany. We prepared 18.4 L of solution at 1.5 wt % gelatin and 15 wt % salt in two steps. First, we dissolved 276 g of pigskin gelatin (280 bloom and 16 mesh, Italgel) in 9 L of water at 60 °C, stirring continuously until the gelatin reached a temperature of 50 °C. Then, in a separate pot, we poured 2760 g of salt into 9.4 L of cold water. To ensure a uniform salt blend in the water, stirring was maintained while the pot with the salty solution was heated to 50 °C. Then, after the two solutions were combined and cooled to approximately 30 °C, we poured the liquid salty-gelatin into a Plexiglas tank with dimensions 40 x 20 cm, reaching a height of 23 cm. The tank was placed in a fridge for 40 hr to solidify the gelatin and cool it down to a stable and homogeneous temperature of 8 °C. We strictly followed the same procedure for preparing all the gelatin blocks used in our experiments to ensure the repeatability of the gelatin characteristics. In addition, the room temperature was maintained below 15 °C to preserve the gelatin's elastic properties during the experiments, which could last up to approximately 3 hr.

Before starting each experiment, we measured the gelatin density and rigidity. Density measurements were performed by weighing and measuring the volume of the gelatin block, obtaining consistent values across all the experiments, ranging between 1101–1129 $\text{kg m}^{-3}$ (Tables 1 and 2). For measuring the rigidity, we followed the procedure described by Pansino & Taisne (2019), using the photoelastic properties of gelatin to measure the velocity of the shear waves. We set two polarized sheets at the front and back of the tank and illuminated the tank to make visible the shear waves produced

by striking the gelatin surface with a plastic spoon. The propagation velocity $v_s$ of the shear waves allowed us to compute the shear modulus $G$:

$$G = \rho v_s^2 \quad (1)$$

considering a Poisson's ratio equal to 0.5 for gelatin (Kavanagh et al. 2013). The values for $G$ are reported in Tables 1 and 2.

We formed fluid-filled cracks of various sizes by injecting different fluid volumes through holes at the bottom of the tank, using a syringe and a metal needle. The needle had a sharp bevelled tip that helped to control the initial crack's orientation and was connected through a silicon hose to the syringe. For each gelatin tank, we performed up to three fluid-filled crack propagation experiments of air and oil injections. Each injection was separated by at least 5 cm from the others to prevent mechanical interactions between a propagating fluid-filled crack and the cuts left behind by previous injections (Le Corvec et al. 2013). Additionally, in two experiments (number 04 and 11), we tested the simultaneous propagation of two oil-filled cracks (cf. Fig. 1a) and observed no differences in the outcomes compared to injections performed in a new—intact—gelatin block or compared to the injections performed successively in the same box. To record the experiments, we set a camera in front of the gelatin tank and used a backlight to enhance the visibility of the fluid-filled cracks. The volume of oil injections ranged from 10 to 50 mL. Given the high viscosity of the silicon oil, a single oil-filled crack propagation experiment could last up to about 3 hr. Thus, we set the camera shoot with a time lapse of 10 s. In contrast, air injections were carried out in volumes ranging from 1 to 10 mL. Due to their low viscosity, these injections display higher velocities, allowing us to record them using videos (Fig. 1c). We ran several experiments injecting the same amount of fluid in order to check the repeatability of the experiments and the stability of our measurements.

### 2.1.2. Measurements of fluid-filled crack aspects and velocity and associated errors

In order to measure the propagation velocity of an ascending fluid-filled crack, we monitored the crack-tip position on its trajectory to the surface. We created videos from image stacks (for oil-filled cracks) and processed them using the free software Tracker (Brown 2012). We tracked the dyke tip position over time manually, and the software provided the velocity between two consecutive positions (Fig. 1c). Additionally, we computed the mean propagation velocity ($v_{\text{mean}}$), its minimum ($v_{\text{min}}$) and maximum values ($v_{\text{max}}$), excluding an initial 'crack formation' phase (marked in grey in Fig. 1c) and the final acceleration when approaching the free surface (Rivalta & Dahm 2006). For oil-filled cracks and air injections up to 5 mL, this phase is characterized by an initially high velocity that quickly decelerates to a more stable trend. In fact, in these cases, the initial crack growth phase occurs at a faster rate than the buoyancy-driven propagation phase. Conversely, the larger air-filled cracks (10 mL) initially display an accelerating velocity. This is due to the fact that these cracks start the buoyant propagation before the end of the injection. As their volume grows, their buoyancy progressively increases, accelerating their propagation velocity until the final volume is reached. In both cases, during these initial phases of deceleration or acceleration, the velocity is significantly influenced by the crack nucleation phase and the rate of fluid injection. These factors are not accounted for in our numerical simulations; hence, their impact cannot be analysed. However, these observations are consistent
<table><caption>**Table 1.** Gelatin properties and experimental results for propagating cracks filled with viscous fluid (silicon oil). The experiments are listed in order of injected volumes (*V*). The first two digits of the ‘Exp. ID’ indicate the ID of the gelatin block, followed by the progressive number of the oil (‘O’) injection, with the exception of the last three rows (bold), listing the averaged values needed to constrain the numerical simulations for 10, 30 and 50 mL, respectively (overlined symbols). Gelatin density (*ρ*), shear modulus (*G*), starting depth (*Z*), as well as physical and geometric crack parameters such as mean-, minimum- and maximum-velocity (*v*<sub>mean</sub>, *v*<sub>min</sub> and *v*<sub>max</sub>), crack opening (*O*), width (*W*) and length (*L*), along with the stress intensity factor estimated from eq. (2) (*K*<sub>j</sub><sup>(3D)</sup>), are reported.</caption><thead><tr><th><i>V</i><br>[mL]</th><th>Exp.<br>ID</th><th><i>ρ</i><br>[kg m<sup>-3</sup>]</th><th><i>G</i><br>[Pa]</th><th><i>Z</i><sub>s</sub><br>[cm]</th><th><i>v</i><sub>mean</sub><br>[mm s<sup>-1</sup>]</th><th><i>v</i><sub>min</sub><br>[mm s<sup>-1</sup>]</th><th><i>v</i><sub>max</sub><br>[mm s<sup>-1</sup>]</th><th><i>O</i><br>[cm]</th><th><i>W</i><br>[cm]</th><th><i>L</i><br>[cm]</th><th><i>K</i><sub>j</sub><sup>(3D)</sup><br>[Pa m<sup>1/2</sup>]</th></tr></thead><tbody><tr><td>10</td><td>04-01</td><td>[11.29 ± 0.34]×10<sup>2</sup></td><td>[1.29 ± 0.13]×10<sup>2</sup></td><td>14.7 ± 0.3</td><td>[1.7 ± 0.3]×10<sup>-2</sup></td><td>[1.0 ± 0.3]×10<sup>-2</sup></td><td>[2.5 ± 0.3]×10<sup>-2</sup></td><td>0.53 ± 0.04</td><td>3.56 ± 0.03</td><td>5.8 ± 0.3</td><td>13 ± 3</td></tr><tr><td>10</td><td>11-01</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>15.6 ± 0.2</td><td>[1.7 ± 0.4]×10<sup>-2</sup></td><td>[0.7 ± 0.4]×10<sup>-2</sup></td><td>[2.6 ± 0.4]×10<sup>-2</sup></td><td>0.49 ± 0.07</td><td>3.38 ± 0.07</td><td>4.9 ± 0.2</td><td>9 ± 2</td></tr><tr><td>10</td><td>11-02</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>15.8 ± 0.2</td><td>[1.7 ± 0.3]×10<sup>-2</sup></td><td>[1.0 ± 0.3]×10<sup>-2</sup></td><td>[2.7 ± 0.3]×10<sup>-2</sup></td><td>0.46 ± 0.03</td><td>3.44 ± 0.12</td><td>5.0 ± 0.2</td><td>9 ± 2</td></tr><tr><td>10</td><td>13-01</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>16.9 ± 0.2</td><td>[1.6 ± 0.4]×10<sup>-2</sup></td><td>[0.5 ± 0.4]×10<sup>-2</sup></td><td>[2.8 ± 0.4]×10<sup>-2</sup></td><td>0.41 ± 0.04</td><td>3.27 ± 0.17</td><td>5.5 ± 0.2</td><td>10 ± 2</td></tr><tr><td>10</td><td>14-01</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.47 ± 0.22]×10<sup>2</sup></td><td>16.2 ± 0.2</td><td>[1.3 ± 0.3]×10<sup>-2</sup></td><td>[0.8 ± 0.3]×10<sup>-2</sup></td><td>[2.1 ± 0.3]×10<sup>-2</sup></td><td>0.37 ± 0.02</td><td>3.31 ± 0.03</td><td>4.9 ± 0.2</td><td>9 ± 2</td></tr><tr><td>20</td><td>05-01</td><td>[11.10 ± 0.29]×10<sup>2</sup></td><td>[1.27 ± 0.11]×10<sup>2</sup></td><td>17.0 ± 0.2</td><td>[4 ± 1]×10<sup>-2</sup></td><td>[2 ± 1]×10<sup>-2</sup></td><td>[7 ± 1]×10<sup>-2</sup></td><td>0.66 ± 0.05</td><td>4.82 ± 0.24</td><td>7.0 ± 0.2</td><td>16 ± 4</td></tr><tr><td>20</td><td>09-02</td><td>[11.01 ± 0.28]×10<sup>2</sup></td><td>[1.34 ± 0.36]×10<sup>2</sup></td><td>15.8 ± 0.2</td><td>[3.8 ± 0.9]×10<sup>-2</sup></td><td>[2.2 ± 0.9]×10<sup>-2</sup></td><td>[5.9 ± 0.9]×10<sup>-2</sup></td><td>0.50 ± 0.05</td><td>4.36 ± 0.28</td><td>6.8 ± 0.2</td><td>13 ± 3</td></tr><tr><td>25</td><td>04-02</td><td>[11.29 ± 0.34]×10<sup>2</sup></td><td>[1.29 ± 0.13]×10<sup>2</sup></td><td>17.9 ± 0.3</td><td>[5 ± 1]×10<sup>-2</sup></td><td>[3 ± 1]×10<sup>-2</sup></td><td>[8 ± 1]×10<sup>-2</sup></td><td>0.79 ± 0.09</td><td>5.60 ± 0.03</td><td>8.3 ± 0.3</td><td>23 ± 5</td></tr><tr><td>30</td><td>08-01</td><td>[11.09 ± 0.29]×10<sup>2</sup></td><td>[1.29 ± 0.20]×10<sup>2</sup></td><td>12.2 ± 0.2</td><td>[5.4 ± 0.9]×10<sup>-2</sup></td><td>[3.1 ± 0.9]×10<sup>-2</sup></td><td>[7.3 ± 0.9]×10<sup>-2</sup></td><td>0.83 ± 0.05</td><td>5.70 ± 0.20</td><td>11.6 ± 0.2</td><td>28 ± 6</td></tr><tr><td>30</td><td>12-01</td><td>[11.06 ± 0.24]×10<sup>2</sup></td><td>[1.29 ± 0.22]×10<sup>2</sup></td><td>13.2 ± 0.1</td><td>[5.2 ± 0.9]×10<sup>-2</sup></td><td>[3.2 ± 0.9]×10<sup>-2</sup></td><td>[7.2 ± 0.9]×10<sup>-2</sup></td><td>0.73 ± 0.02</td><td>5.91 ± 0.30</td><td>9.0 ± 0.1</td><td>22 ± 4</td></tr><tr><td>30</td><td>12-02</td><td>[11.06 ± 0.24]×10<sup>2</sup></td><td>[1.29 ± 0.22]×10<sup>2</sup></td><td>11.7 ± 0.1</td><td>[5.4 ± 0.9]×10<sup>-2</sup></td><td>[3.2 ± 0.9]×10<sup>-2</sup></td><td>[7.3 ± 0.9]×10<sup>-2</sup></td><td>0.73 ± 0.02</td><td>5.42 ± 0.06</td><td>8.5 ± 0.1</td><td>20 ± 4</td></tr><tr><td>30</td><td>13-02</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>19.1 ± 0.2</td><td>[5.6 ± 0.7]×10<sup>-2</sup></td><td>[3.1 ± 0.7]×10<sup>-2</sup></td><td>[8.1 ± 0.7]×10<sup>-2</sup></td><td>0.71 ± 0.05</td><td>5.21 ± 0.54</td><td>8.6 ± 0.2</td><td>19 ± 4</td></tr><tr><td>30</td><td>14-02</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.47 ± 0.22]×10<sup>2</sup></td><td>15.6 ± 0.2</td><td>[4.8 ± 0.9]×10<sup>-2</sup></td><td>[2.9 ± 0.9]×10<sup>-2</sup></td><td>[7.4 ± 0.9]×10<sup>-2</sup></td><td>0.88 ± 0.01</td><td>5.74 ± 0.03</td><td>7.0 ± 0.2</td><td>17 ± 4</td></tr><tr><td>35</td><td>04-03</td><td>[11.29 ± 0.34]×10<sup>2</sup></td><td>[1.29 ± 0.13]×10<sup>2</sup></td><td>18.8 ± 0.3</td><td>[6 ± 1]×10<sup>-2</sup></td><td>[3 ± 1]×10<sup>-2</sup></td><td>[8 ± 1]×10<sup>-2</sup></td><td>0.98 ± 0.02</td><td>5.89 ± 0.01</td><td>10.6 ± 0.3</td><td>30 ± 7</td></tr><tr><td>35</td><td>09-03</td><td>[11.01 ± 0.28]×10<sup>2</sup></td><td>[1.34 ± 0.36]×10<sup>2</sup></td><td>19.2 ± 0.2</td><td>[6.1 ± 0.8]×10<sup>-2</sup></td><td>[4.4 ± 0.8]×10<sup>-2</sup></td><td>[7.8 ± 0.8]×10<sup>-2</sup></td><td>0.81 ± 0.04</td><td>5.83 ± 0.06</td><td>11.8 ± 0.2</td><td>27 ± 6</td></tr><tr><td>40</td><td>05-02</td><td>[11.10 ± 0.29]×10<sup>2</sup></td><td>[1.27 ± 0.11]×10<sup>2</sup></td><td>16.3 ± 0.2</td><td>[6.8 ± 0.9]×10<sup>-2</sup></td><td>[4.1 ± 0.9]×10<sup>-2</sup></td><td>[9.4 ± 0.9]×10<sup>-2</sup></td><td>0.91 ± 0.09</td><td>6.37 ± 0.10</td><td>13.0 ± 0.2</td><td>33 ± 7</td></tr><tr><td>50</td><td>08-02</td><td>[11.09 ± 0.29]×10<sup>2</sup></td><td>[1.29 ± 0.20]×10<sup>2</sup></td><td>15.9 ± 0.2</td><td>[8 ± 2]×10<sup>-2</sup></td><td>[4 ± 2]×10<sup>-2</sup></td><td>[1.3 ± 0.2]×10<sup>-1</sup></td><td>0.97 ± 0.03</td><td>6.78 ± 0.29</td><td>14.1 ± 0.2</td><td>37 ± 8</td></tr><tr><td>50</td><td>09-01</td><td>[11.01 ± 0.28]×10<sup>2</sup></td><td>[1.34 ± 0.36]×10<sup>2</sup></td><td>17.7 ± 0.2</td><td>[8 ± 2]×10<sup>-2</sup></td><td>[4 ± 2]×10<sup>-2</sup></td><td>[1.0 ± 0.2]×10<sup>-1</sup></td><td>0.98 ± 0.07</td><td>5.92 ± 0.19</td><td>18.0 ± 0.2</td><td>42 ± 10</td></tr><tr><td>50</td><td>10-01</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>20.4 ± 0.2</td><td>[8 ± 1]×10<sup>-2</sup></td><td>[4 ± 2]×10<sup>-2</sup></td><td>[1.1 ± 0.2]×10<sup>-1</sup></td><td>1.00 ± 0.06</td><td>6.92 ± 0.23</td><td>15.4 ± 0.2</td><td>40 ± 9</td></tr><tr><td>50</td><td>13-03</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>15.7 ± 0.2</td><td>[10 ± 1]×10<sup>-2</sup></td><td>[8 ± 1]×10<sup>-2</sup></td><td>[1.2 ± 0.1]×10<sup>-1</sup></td><td>1.00 ± 0.04</td><td>5.76 ± 0.35</td><td>15.9 ± 0.2</td><td>38 ± 8</td></tr><tr><td>50</td><td>14-03</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.47 ± 0.22]×10<sup>2</sup></td><td>12.8 ± 0.2</td><td>[7 ± 2]×10<sup>-2</sup></td><td>[4 ± 2]×10<sup>-2</sup></td><td>[1.3 ± 0.2]×10<sup>-1</sup></td><td>0.93 ± 0.02</td><td>7.58 ± 0.12</td><td>10.3 ± 0.2</td><td>28 ± 6</td></tr><tr><td><i>V</i><br>[mL]</td><td>Sim.<br>ID</td><td><i>ρ</i><br>[kg m<sup>-3</sup>]</td><td><i>G̃</i><br>[Pa]</td><td>—</td><td><i>v̇</i><sub>mean</sub><br>[mm s<sup>-1</sup>]</td><td><i>v̇</i><sub>min</sub><br>[mm s<sup>-1</sup>]</td><td><i>v̇</i><sub>max</sub><br>[mm s<sup>-1</sup>]</td><td><i>Ȯ</i><br>[cm]</td><td>—</td><td>—</td><td>—</td></tr><tr><td>10</td><td>10-AV</td><td>[11.10 ± 0.29]×10<sup>2</sup></td><td>[1.35 ± 0.22]×10<sup>2</sup></td><td>—</td><td>1.6 × 10<sup>-2</sup></td><td>0.8 × 10<sup>-2</sup></td><td>2.5 × 10<sup>-2</sup></td><td>0.45 ± 0.06</td><td>—</td><td>—</td><td>—</td></tr><tr><td>30</td><td>30-AV</td><td>[11.07 ± 0.27]×10<sup>2</sup></td><td>[1.33 ± 0.22]×10<sup>2</sup></td><td>—</td><td>5.3 × 10<sup>-2</sup></td><td>3.1 × 10<sup>-2</sup></td><td>7.5 × 10<sup>-2</sup></td><td>0.78 ± 0.07</td><td>—</td><td>—</td><td>—</td></tr><tr><td>50</td><td>50-AV</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.34 ± 0.24]×10<sup>2</sup></td><td>—</td><td>8 × 10<sup>-2</sup></td><td>5 × 10<sup>-2</sup></td><td>1.2 × 10<sup>-1</sup></td><td>0.98 ± 0.05</td><td>—</td><td>—</td><td>—</td></tr></tbody></table>
<table><caption>**Table 2.** Gelatin properties and experimental results for propagating air-filled cracks. The experiments are listed in order of injected volumes (*V*). The first two digits of the ‘Exp. ID’ indicate the ID of the gelatin block, followed by the progressive number of the air (‘A’) injection. The gelatin density (*ρ*), shear modulus (*G*), starting depth (*Z*<sub>s</sub>), as well as physical and geometric crack parameters such as mean-, minimum- and maximum-velocity (*v*<sub>mean</sub>, *v*<sub>min</sub> and *v*<sub>max</sub>), crack opening (*O*), width (*W*) and length (*L*), along with the stress intensity factor estimated from eq. (2) (*K*<sub>f</sub><sup>(3D)</sup>), are reported.</caption><thead><tr><th><i>V</i><br>[mL]</th><th>Exp.<br>ID</th><th><i>ρ</i><br>[kg m<sup>-3</sup>]</th><th><i>G</i><br>[Pa]</th><th><i>Z</i><sub>s</sub><br>[cm]</th><th><i>v</i><sub>mean</sub><br>[mm s<sup>-1</sup>]</th><th><i>v</i><sub>min</sub><br>[mm s<sup>-1</sup>]</th><th><i>v</i><sub>max</sub><br>[mm s<sup>-1</sup>]</th><th><i>O</i><br>[cm]</th><th><i>W</i><br>[cm]</th><th><i>L</i><br>[cm]</th><th><i>K</i><sub>f</sub><sup>(3D)</sup><br>[Pa·m<sup>1/2</sup>]</th></tr></thead><tbody><tr><td>1</td><td>10-A1</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>9.7 ± 0.2</td><td>2.29 ± 0.15</td><td>2.25 ± 0.15</td><td>2.38 ± 0.15</td><td>0.65 ± 0.02</td><td>1.51 ± 0.06</td><td>2.3 ± 0.2</td><td>25 ± 2</td></tr><tr><td>1</td><td>11-A1</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>9.3 ± 0.2</td><td>1.98 ± 0.17</td><td>1.89 ± 0.17</td><td>2.04 ± 0.17</td><td>0.66 ± 0.06</td><td>1.40 ± 0.08</td><td>2.4 ± 0.2</td><td>25 ± 2</td></tr><tr><td>1</td><td>12-A1</td><td>[11.06 ± 0.24]×10<sup>2</sup></td><td>[1.29 ± 0.22]×10<sup>2</sup></td><td>17.1 ± 0.1</td><td>1.7 ± 0.4</td><td>0.9 ± 0.4</td><td>2.1 ± 0.4</td><td>0.62 ± 0.02</td><td>1.41 ± 0.06</td><td>2.4 ± 0.1</td><td>25 ± 1</td></tr><tr><td>2</td><td>11-A2</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>19.9 ± 0.2</td><td>6 ± 1</td><td>5 ± 1</td><td>8 ± 1</td><td>0.75 ± 0.02</td><td>1.77 ± 0.05</td><td>2.5 ± 0.2</td><td>29 ± 2</td></tr><tr><td>2</td><td>12-A2</td><td>[11.06 ± 0.24]×10<sup>2</sup></td><td>[1.29 ± 0.22]×10<sup>2</sup></td><td>20.3 ± 0.1</td><td>7 ± 1</td><td>5 ± 1</td><td>8 ± 1</td><td>0.81 ± 0.03</td><td>1.91 ± 0.07</td><td>2.7 ± 0.1</td><td>32 ± 2</td></tr><tr><td>2</td><td>13-A1</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>20.4 ± 0.2</td><td>7 ± 1</td><td>5 ± 1</td><td>9 ± 1</td><td>0.97 ± 0.02</td><td>1.91 ± 0.03</td><td>2.8 ± 0.2</td><td>34 ± 3</td></tr><tr><td>3</td><td>10-A2</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>19.2 ± 0.2</td><td>12 ± 1</td><td>9 ± 1</td><td>14 ± 1</td><td>1.05 ± 0.03</td><td>2.27 ± 0.02</td><td>2.8 ± 0.2</td><td>37 ± 3</td></tr><tr><td>3</td><td>11-A3</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>19.2 ± 0.2</td><td>11 ± 2</td><td>9 ± 2</td><td>14 ± 2</td><td>1.01 ± 0.02</td><td>2.25 ± 0.02</td><td>2.9 ± 0.2</td><td>38 ± 3</td></tr><tr><td>3</td><td>12-A3</td><td>[11.06 ± 0.24]×10<sup>2</sup></td><td>[1.29 ± 0.22]×10<sup>2</sup></td><td>20.3 ± 0.1</td><td>11 ± 1</td><td>9 ± 1</td><td>13 ± 1</td><td>1.10 ± 0.03</td><td>2.31 ± 0.06</td><td>2.9 ± 0.1</td><td>38 ± 2</td></tr><tr><td>4</td><td>13-A2</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>19.3 ± 0.2</td><td>15 ± 2</td><td>13 ± 2</td><td>18 ± 2</td><td>1.35 ± 0.06</td><td>2.46 ± 0.01</td><td>3.1 ± 0.2</td><td>42 ± 3</td></tr><tr><td>5</td><td>10-A3</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>18.4 ± 0.2</td><td>19 ± 2</td><td>17 ± 2</td><td>22 ± 2</td><td>1.35 ± 0.06</td><td>2.76 ± 0.07</td><td>3.4 ± 0.2</td><td>49 ± 3</td></tr><tr><td>5</td><td>11-A4</td><td>[11.05 ± 0.28]×10<sup>2</sup></td><td>[1.32 ± 0.27]×10<sup>2</sup></td><td>18.3 ± 0.2</td><td>18 ± 2</td><td>15 ± 2</td><td>22 ± 2</td><td>1.32 ± 0.01</td><td>2.77 ± 0.03</td><td>3.4 ± 0.2</td><td>49 ± 3</td></tr><tr><td>5</td><td>13-A3</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.33 ± 0.21]×10<sup>2</sup></td><td>19.0 ± 0.2</td><td>19 ± 2</td><td>16 ± 2</td><td>22 ± 2</td><td>1.39 ± 0.02</td><td>2.81 ± 0.04</td><td>3.3 ± 0.2</td><td>48 ± 3</td></tr><tr><td>10</td><td>10-A4</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>15.4 ± 0.2</td><td>36 ± 2</td><td>32 ± 2</td><td>39 ± 2</td><td>1.79 ± 0.09</td><td>3.57 ± 0.03</td><td>3.9 ± 0.2</td><td>64 ± 4</td></tr><tr><td>10</td><td>10-A5</td><td>[11.06 ± 0.28]×10<sup>2</sup></td><td>[1.25 ± 0.23]×10<sup>2</sup></td><td>14.9 ± 0.2</td><td>36 ± 2</td><td>34 ± 2</td><td>37 ± 2</td><td>1.69 ± 0.08</td><td>3.44 ± 0.07</td><td>3.9 ± 0.2</td><td>63 ± 4</td></tr></tbody></table>

<figure><img src="image_2.png" /><figcaption>**Figure 1.** (a) Laboratory setup and image of two parallel oil-filled crack propagation corresponding to experiments 11-O2 (left) and 11-O1 (right) (see Table 1). (b) Top-view image of experiment 11-O1 highlighting the crack width (*W*) and opening (*O*) just before the arrival at surface. (c) Crack propagation velocity (*v*) against depth (*Z*) for a 10 mL oil-filled crack experiment (11-O1) and a 10 mL air-filled crack (10-A5). Note that the velocity range for an oil-filled crack is 3 orders of magnitude smaller than for an air-filled crack for the same injected volume. The opposite behaviours (deceleration and acceleration) observed for the 10 mL oil- and air-filled cracks during the ‘initial phase’ (marked in grey) are mainly due to the large buoyancy difference of these fluids and are discussed in Section 2.1.2.</figcaption></figure>

with the theoretical considerations and the numerical results obtained by Möri & Lecampion (2022, 2023), who highlight similar dynamics in the transition from initial radial to buoyancy-driven propagation. For each experiment, we identified the depth reached by the dyke tip once the crack formation phase was over (*Z*<sub>s</sub> in Table 1). Crack characteristics, including its maximum opening (*O*), head width (*W*) and length (*L*), were measured from the recorded images using Tracker. We used images from the front camera to measure *L*, while *O* and *W* were measured from images taken from a top view. To mitigate the influence of light refraction in the gelatin block (Furst *et al.* 2024), measurements from top-view images were

taken when the crack reached about 2 cm from the surface. For consistency, the same was applied also to length measurements (cf. Fig. 1b). Images were calibrated using transparent graduated sheets placed on the gelatin surface and the tank walls. Some cracks formed with a slight dip angle with respect to the vertical (e.g. 11-O1 in Fig. 1), which was due to the orientation of the injection needle (controlling the initial fracture orientation). Consistently producing exactly vertical cracks, as well as cracks which are perpendicular to the front camera, remains challenging. Therefore, we sometimes needed to change the view angle of the camera to capture the crack opening and width. This procedure introduced a distortion that had
to be taken into account during the image calibration in order to measure O and W.

Uncertainties (Tables 1 and 2) were estimated from direct measurement errors, taking into account the sensitivity of our instruments, but also other possible error sources. In particular, the errors associated with measurements taken on recorded images strongly depend on the quality of the image, that is, its sharpness, resolution and possible distortion effects. For instance, when measuring the length L of a crack on a digital image, we considered its maximum and minimum estimates in pixels, accounting for the thickness of the line marking the crack borders on our recordings. The associated error was taken as the difference between these two estimates. As a consequence, the sharper the recorded image, the thinner the line marking the crack border, and the smaller the error. Additionally, when scaling the pixels of an image, we accounted for distortion effects at its edges using transparent graduated sheets placed on the sides of the gelatin tank. Similarly, for the manual pick of the crack front, we considered the midpoint of the crack borderline. These picks were used to compute the crack propagation velocity using the software ‘Tracker’. The errors associated with the velocity were computed as their standard deviation within a time interval with relatively stable average velocity. Eventually, the errors associated with derived quantities were computed by combining the relative errors associated with direct measurements (for instance, the error associated with the gelatin density in Tables 1 and 2, was computed using $\sigma_\rho/\rho = \sigma_M/M + \sigma_\nu/V$).

### 2.1.3. Fracture toughness computation

The fracture toughness $K_c$ is a parameter used to characterize the resistance to the fracturing of a brittle material. A fluid-filled crack is stable (does not grow larger) when the stress intensity factor at the crack tip is lower than the fracture toughness of the host solid. The stress intensity factor $K_I$ depends on the fluid excess pressure with respect to the confining stress and the shape of the crack. 3-D and 2-D stress intensity factor estimates can be provided by the following static crack formulae (Secor & Pollard 1975; Heimpel & Olson 1994, respectively):

$$K_I^{(3D)} = \sqrt{\frac{2}{\pi}} \Delta\rho \cdot g \cdot L \sqrt{W}, \qquad (2)$$

$$K_I^{(2D)} = \sqrt{\frac{\pi}{8}} \Delta\rho \cdot g \cdot L^{3/2} \qquad (3)$$

where $\Delta\rho$ is the density contrast between the fluid and the solid, and g is the acceleration due to gravity. Using eq. (2), we computed $K_I^{(3D)}$ for fluid-filled cracks with different volumes, propagating at different velocities, and provided an estimate of the gelatin fracture toughness $K_c$ by extrapolating $K_I^{(3D)}(v)$ for the velocity approaching zero (Smittarello et al. 2021). In fact, even though $K_I^{(3D)}$ does not explicitly depend on the fluid-filled crack propagation velocity (eq. 2), it depends on the injected fluid volume (affecting L and W). In our experiments, a larger fluid volume implies a larger overpressure and an increased crack propagation velocity. This is generally true for toughness-dominated buoyant cracks, while during the transition between viscosity- to toughness-dominated regimes, the velocity increases even though the crack head volume actually decreases (Möri & Lecampion 2022). However, within the range of our experiments, we do not observe such a transition in the propagation regimes.

In the assumption that $K_I = K_{eff}(v)$ during the fluid-filled crack propagation, extrapolating $K_I$ values for $v \to 0$ is equivalent to

computing $K_c = K_{eff}(v \to 0)$. Note that more general frameworks for the dynamics of fracture propagation consider different (more complex) velocity dependencies for $K_I(v)$ and $K_{eff}(v)$ (Bhat et al. 2012), as compared to the relations used for quasi-static fluid-filled fractures.

Additionally, eq. (2) applies to stationary cracks, and it can be used for quasi-static moving cracks only when viscous forces do not alter their shape. This is true for air-filled cracks, as their shape is indistinguishable from the static one (Dahm 2000a). Several previous studies that used air- and oil-filled crack injections in gelatin (e.g. Takada 1990; Heimpel & Olson 1994; Smittarello et al. 2021) neglected the effect of fluid viscosity, in the assumption that either the fluid viscosity was sufficiently low (for air-filled cracks), or the fluid-filled crack propagation velocity was slow enough (for oil-filled cracks).

Here, we will first compute $K_I^{(3D)}$ for air- and oil-filled crack propagation experiments using eq. (2), and we will estimate $K_c$ for both sets of experiments separately. However, we will additionally provide independent estimates of $K_I^{(3D)}$ for the oil-filled cracks, accounting for the oil viscosity, since the shape of these cracks in our experiments differs from the static case. In order to do so, we will make use of 2-D numerical simulations applied to the oil-filled cracks. Finally, to compare results from 3-D experiments and 2-D simulations, we introduce a method to convert 3-D and 2-D stress intensity factor estimates.

Smittarello et al. (2021) used 3-D and 2-D stress intensity factor formulae to derive a relation between them, combining eqs (2) and (3):

$$K_I^{(3D)} = \frac{4}{\pi} \sqrt{\frac{W}{L}} \cdot K_I^{(2D)}. \qquad (4)$$

This relation will be used to consistently show on the same plots results from (3-D) laboratory experiments and (2-D) numerical simulations. In particular, we will use W/L average values for 10, 30 and 50 mL oil injections when converting to 3-D, the 2-D estimates obtained with the numerical simulations of those experiments, as the numerical model is in plane strain and does not provide W values. In the Discussion (Section 4.2), we will provide further details about W/L measurements and their use in eq. (4) to convert between 2-D and 3-D estimates.

## 2.2. Numerical model for viscous fluids

In order to account for the effect of oil viscosity on the stress intensity factor calculation and fracture toughness estimates, we used the numerical model developed by Furst et al. (2023). We did not use the numerical model to simulate the air-filled crack propagation experiments as the viscosity of air is so low that the dynamic-crack shape and pressure profile coincides with the analytical static-crack solution, confirming the validity of eq. (2) for our air-filled cracks. The 2-D model developed by Furst et al. (2023) integrates fracture mechanics with fluid flow equations to simulate the evolution of a fluid-filled crack in terms of its shape and propagation velocity. This model refers to a vertical cross-section of the fracture, discretized using dislocation elements in plane strain approximation (e.g. the fracture is considered to extend infinitely in the out-of-plane direction). The model simulates a fracture filled with a compressible fluid flowing between the crack walls in a Hagen-Poiseuille flow regime. The crack propagation velocity (v) is determined through an energy balance criterion, performed over the entire crack length, and considering different energy contributions: potential energy variations ($\Delta E$), which include elastic and gravitational energy changes
associated with the upward propagation of a buoyant fluid-filled crack; fracture energy ($\Delta E_f$), the energy per unit surface required to fracture the host solid at the crack's tip, multiplied by the fracture elongation ($l$); and viscous energy dissipation ($\Delta E_v$), arising from the flow of viscous fluid between the crack walls. Due to the plane strain (2-D) approximation, all energy terms are expressed as energy per unit length, also referred to as 'energy density' (Rice 1968), and the energy budget equation is:

$$ \Delta E = \Delta E_f \cdot l + \Delta E_v \quad (5) $$

with $\Delta E_v = (v \cdot \eta)D'$, where $\eta$ is the fluid viscosity and $D'$ is a function that depends on the crack geometry and magma flow (Furst et al. 2023, where $D = \eta \cdot D'$). The energy budget equation (5) introduces a direct constraint on the fluid-filled crack propagation velocity:

$$ v = \frac{\Delta E - \Delta E_f \cdot l}{\eta \cdot D'} \quad (6) $$

Note that $\Delta E$ and $D'$ intrinsically depend on $v$, and eq. (6) must be solved iteratively (Furst et al. 2023).

$\Delta E_f$ depends on the fracture toughness $K_c$ of an elastic material through the following relation (e.g. Rice 1968; Dahm 2000b):

$$ \Delta E_f = K_c^2 \frac{1 - nu}{2G} \quad (7) $$

In the numerical approach (Furst et al. 2023), the fracture energy per unit surface $\Delta E_f$ is constant during the crack growth; therefore, $K_c$ also is assumed to remain constant. As a consequence, for a given value of $K_c$ such that $\Delta E_f \cdot l < \Delta E$, the viscous flow within the crack represents the only limiting factor for the propagation velocity of the fluid-filled crack. However, if the viscosity approaches zero, eq. (6) is not suitable for computing $v$, as the velocity would grow to unrealistic values (as in the case of air-filled cracks), making this model inapplicable for the computation of the propagation velocity of low viscosity (fracture-dominated) fluid-filled cracks (Furst et al. 2023). Therefore, here, we performed numerical runs directly constrained by the oil-filled crack experiments.

The model requires several inputs: the mechanical properties of the fluid and the embedding medium ($\rho$, $G$, $nu$, $\rho_f$ and $\eta$), an initial velocity guess ($v_i$, which is a numerical parameter needed to reach the algorithm convergence, without affecting the final results), the cross-sectional area of the fluid-filled crack at a reference pressure ($A_0$), and the fracture energy ($\Delta E_f$). As output, the simulations provide the crack shape and its velocity profile during propagation. Although we could obtain precise measurements for most input parameters, such as $\rho$, $G$, $nu$, $\rho_f$, $\eta$ and the velocity guess (taken as the average observed velocity $\bar{v}_{mean}$, Table 1), two parameters remained unknown or poorly constrained: the fracture energy ($\Delta E_f$) and the cross-sectional area of the crack ($A_0$), respectively. In fact, measuring the crack area is particularly challenging as the cracks are thin, and their orientation is not always perpendicular to the camera view (Smittarello et al. 2021). Also, estimating the 2-D cross-section using the actual 3-D injected volume is not straightforward, as we could only measure the crack head length $L$, its maximum opening $O$, and width $W$, but not their distributed values over the crack surface. As a consequence, following Furst et al. (2024), for each simulation we looked for the values of $A_0$ and $\Delta E_f$ that provided the best simultaneous fit for both the maximum crack opening ($O$) and the observed propagation velocity (considering the $v_{mean}$, $v_{min}$ and $v_{max}$ values averaged over all the experiments with the same volumes indicated as overlined variables, last three rows of Table 1). In fact, this is the combination of parameters that are deemed

the most reliable quantities to constrain the Furst et al. (2023) model, as it has been shown for similar experiments in Furst et al. (2024).

In order to perform the fit between the simulated and observed velocity, here we followed a simple trial-and-error approach, varying $\Delta E_f$ and visually adjusting the simulated velocity profiles—excluding the initial rapid velocity decay—to the measured velocities. This procedure worked efficiently as different $\Delta E_f$ values would essentially shift the simulated velocity profile towards higher or lower velocity values.

Given $\Delta E_f$, eq. (7) would then provide estimates of the gelatin fracture toughness $K_c$ for each simulated experiment. If the fluid-filled fracture propagation occurs as a purely brittle elastic process, $\Delta E_f$ values (and the corresponding estimates of $K_c$) must be equal for all simulations. Additionally, $K_c$ estimates from oil- and air-filled crack experiments should also be equal once we accounted for the oil viscosity and corrected for the 2-D approximation of the numerical simulations (eq. 4, taking $K_c = K_l$).

However, as we will show in the results section, our numerical simulations provide us with different $\Delta E_f$ estimates for different oil-filled crack simulations propagating with different velocities. Because of this, we implemented a modified version of the numerical model, introducing a velocity-dependent relation for the fracture energy $\Delta E_f(v)$. Such relation will be computed by interpolating $\Delta E_f$ values obtained from the different simulations characterized by different propagation velocities. This is equivalent to considering effective fracture toughness values $K_{eff}(v)$, instead of a constant $K_c$, and a fracture propagation condition $K_l = K_{eff}(v)$. Consistently, we should replace eq. (7) with:

$$ \Delta E_f(v) = K_{eff}^2(v) \frac{1 - nu}{2G} \quad (8) $$

In the next section, we will also present results obtained by introducing velocity-dependent fracture energy $\Delta E_f(v)$ within the numerical code for fluid-filled fracture propagation by Furst et al. (2023). This was achieved by implementing an empirical expression for $\Delta E_f(v)$ in the numerical code provided as input. At each time step, the energy balance (eq. 5) is solved iteratively, starting with the velocity value obtained by the crack growth at the previous time step and updating it along with all energy terms until convergence is reached.

Eventually, as our numerical model does not provide direct estimates of $K_l$ values, estimates of $K_l(v)$ for the oil-filled cracks will be provided by the values of $\Delta E_f(v)$ obtained with the numerical simulations, using eq. (8), and the propagation condition $K_l = K_{eff}(v)$. Note that these estimates consider the effect of fluid viscosity, as numerical simulations account for the viscous dissipation in the energy budget (eq. 5).

### 3. RESULTS

In Tables 1 and 2, we present our measurements and results for 21 silicon oil injection experiments, with volumes ranging from 10 to 50 mL and 15 air injections, with volumes ranging from 1 to 10 mL. All injections were performed within eight different gelatin blocks with very similar elastic properties. Oil-filled cracks display stress intensity factors ranging between $9 < K_l^{(oil)} < 41$ Pa·m$^{1/2}$, associated with average velocities in the order of $10^{-2}$ mm s$^{-1}$. For air-filled cracks, we obtain $24 < K_l^{(air)} < 64$ Pa·m$^{1/2}$ for cracks propagating with a velocity between 1 and 36 mm s$^{-1}$. In Fig. 2, we plotted the stress intensity factor against the mean propagation velocity ($v_{mean}$) for each oil injection experiment (Fig. 2a), and for the air injections
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Stress intensity factor $K_I$ plotted against mean velocities ($v_{mean}$) obtained from fluid-filled crack propagating in gelatin. Note that in the assumption of crack propagation with $K_I = K_{eff}(v)$, the y-axis can also be interpreted as the gelatin's effective fracture toughness. (a) Cracks filled with viscous fluid (silicon oil) and volumes ranging from 10 to 50 mL: 10 mL in squares, 30 mL in triangles, 50 mL in diamonds and others in circles (Table 1). (b) Cracks filled with non-viscous fluid (air) and volumes ranging from 1 to 10 mL (Table 2). For each set of experiments, the dashed lines represent the linear regression of $K_I$ as a function of $v_{mean}$: $K_I^{(\text{oil})} = 404.9v + 2.2$ and $K_I^{(\text{air})} = 1.2v + 24.0$.</figcaption></figure>

(Fig. 2b). The dashed lines represent the linear regression of $K_I$ as a function of $v_{mean}$. Assuming that the crack propagation occurs with a stress intensity factor equal to fracture toughness, $K_I = K_{eff}(v)$, the y-intercept provides an estimate of $K_c = K_{eff}(v \to 0)$ for both sets of experiments. Although the characteristics of the gelatin blocks are very similar across all experiments, the estimates of $K_c$ obtained with the two sets of experiments differ by more than one order of magnitude: $K_c^{(\text{oil})} = 2.2 \text{ Pa} \cdot \text{m}^{1/2}$ and $K_c^{(\text{air})} = 24.0 \text{ Pa} \cdot \text{m}^{1/2}$ for oil and air, respectively. This shows that the two linear regressions for $K_I(v)$, obtained over very different ranges of crack propagation velocities (for oil- and air-filled cracks), are not adequate to extrapolate consistent estimates of gelatin strength. Rather, considering $K_I = K_{eff}(v)$, the two regressions displayed in Fig. 2 suggest that the effective fracture toughness should be characterized by different velocity-dependent relations over such largely different velocity ranges.

Additionally, in Fig. 2, we used eq. (2) to compute $K_I$. However, eq. (2) neglects the effect of fluid viscosity on the fluid pressure profile and consequently disregards the viscous energy dissipation due to the internal friction caused by the viscous shear flow of a fluid, which may be important for the estimation of $K_I$ for oil-filled cracks.

In order to overcome this limitation and find improved relations for $K_{eff}(v)$, we used the numerical model from Furst *et al.* (2023), which accounts for the effect of fluid viscous motion, and performed three simulations of oil-filled crack propagation corresponding to the experiments with 10, 30 and 50 mL oil injections. The simulations were constrained by averaging the geometrical and physical parameters measured for 10, 30 and 50 mL oil-injection experiments, listed in Table 1, and labelled as 10-AV, 30-AV and 50-AV. The other input parameters are the oil density $\rho_{oil} = 980 \text{ kg m}^{-3}$, and viscosity $\eta = 9.7 \text{ Pa} \cdot \text{s}$, both constant for all simulations, and the initial crack cross-sections: $A_0 = 2.4, 4.9$ and $6.7 \text{ cm}^2$, respectively. Note that these values of $A_0$ were constrained using the crack maximum openings *O* and velocity *v*, as we explained in Section 2.2. The approximate 3-D volume corresponding to $A_0$ can be estimated by multiplying $A_0$ times the measured crack width *W* (not provided by the 2-D model). Even though this represents quite a rough estimate, using *W* mean values for 10, 30 and 50 mL oil experiments, we obtained crack volumes which are compatible with the actual

injected volumes, with relative differences ranging between 10 and 20% (see Supporting Information, Table S1). In Fig. 3, we display the full velocity profiles (grey lines) as recorded throughout the entire duration of each oil-injection experiment, together with the velocity obtained with the numerical simulations (bold lines), for 10, 30 and 50 mL injected volumes (Figs 3a, b and c, respectively). The observed velocity profiles show that during the initial phase of an oil injection experiment, the velocity of the propagating tip of the oil-filled cracks decreases relatively fast. This phase is governed by the formation of the fluid-filled fracture, and the velocity observed during this initial phase of propagation is also strongly affected by the oil injection rate, as we will further discuss in Section 4. However, our numerical simulations do not account for the fluid injection rate because each simulation starts with all the fluid already in place within an initial crack length *L* (cf. Table 1, rows 10-AV, 30-AV and 50-AV). Therefore, we discard from our analysis the initial part of the velocity profiles below the depth $Z_s$, identified by the fast decrease in the observed propagation velocity (areas of fast velocity decay in Fig. 3).

In Fig. 3, the horizontal lines represent the values of the mean, maximum and minimum velocities averaged over each set of experiments for depths above $Z_s$: $\bar{v}_{mean}$ (solid line), $\bar{v}_{min}$ and $\bar{v}_{max}$ (dashed lines), respectively (Table 1). For each simulation, we looked for the value of the gelatin fracture energy $\Delta E_f$ that provides a good overlap between the simulated and observed velocity profiles. In fact, decreasing $\Delta E_f$ would essentially shift the simulated velocity profile towards higher velocities and vice versa. Following the initial phase of fast velocity decay, both laboratory experiments and numerical simulations show velocity profiles that exhibit a stable decelerating trend until approaching the free surface, where acceleration is observed. A progressive velocity decrease, following the initial crack growth phase, is also predicted by the lubrication theory (Spence & Turcotte 1990). However, our numerical simulations—that converge to (Spence & Turcotte 1990) solution in the viscous-dominated regime (Furst *et al.* 2023)—display a faster velocity decay than the experiments, and a larger final acceleration (cf. Fig. 3). These differences cannot be reduced when using a constant value for $\Delta E_f$, as we will show in more detail later in this section.

In Fig. 3d, we plot the values of $\Delta E_f$ obtained from each of the three simulations as a function of their average velocities, and
<figure><img src="image_4.png" /><figcaption>**Figure 3.** (a), (b) and (c) Velocity variations as a function of depth for 10, 30 and 50 mL oil injections, respectively. Experimental records are plotted in grey lines in the background; bold lines are the numerical simulations; horizontal solid lines and horizontal dashed lines are the mean, maximum and minimum velocities averaged over all experiments with 10, 30 and 50 mL injections ($\bar{v}_{mean}$, $\bar{v}_{min}$ and $\bar{v}_{max}$), (Table 1). The vertical dotted lines represent the end of the initial growing phase of the crack (depth $Z > Z_s$). (d) Fracture energy estimated from the numerical simulations versus average propagation velocity ($\bar{v}_{mean}$). Square, triangle and star markers represent 10, 30 and 50 mL numerical equivalent cracks. The linear relationship depends on the propagation velocity $\bar{v}_{mean}$ expressed in mm s⁻¹.</figcaption></figure>

we display the best-fitting linear trend across these points. Our numerical simulations suggest that a simple linear relation $\Delta E_f(v) = 6.14 v + 0.09$ Pa·m (which corresponds to a square root relation for $K_{eff}(v)$, using eq. 8), may be sufficient to describe fracture energy variations within the velocity range explored in our oil-filled crack experiments (Fig. 3d). In general, this supports the idea that fracture energy (as well as gelatin fracture toughness) is an effective parameter that depends on the fracture velocity.

Therefore, we modified the numerical model for fluid-filled fracture propagation accounting for the linear relationship between $\Delta E_f$ and $v$ displayed in Fig. 3d. In Fig. 4, we show the results from the new numerical simulations performed for the 10, 30 and 50 mL oil injections, with $\Delta E_f(v) = 6.14 v + 0.09$ Pa·m. The fit between model results and observed velocity profiles improves considerably within the range of model applicability (e.g. when the initial fast velocity decay ends). Particularly, the final acceleration due to the free surface obtained with the new simulations reproduces the experimental observations much better.

Since the numerical simulations supported the hypothesis of velocity-dependent fracture energy $\Delta E_f(v)$, we used $\Delta E_f(v = \bar{v}_{mean})$ values and eq. (8), along with the propagation condition $K_I = K_{eff}(v)$, to provide new estimates of the stress intensity factor $K_I$ for the oil-filled cracks. As the numerical simulations are 2-D, these estimates of $K_I$ are also in 2-D. To make these values directly comparable to previous 3-D estimates, we used eq. (4), considering average $W/L$

values for each volume experiment ($W/L = 0.65, 0.64$ and $0.47$ for 10, 30 and 50 mL, respectively), obtaining the empty symbols displayed in Fig. 5a. In Fig. 5a, we also display $K_I^{(oil)}$ values previously computed using eq. (2) (which does not account for the effect of fluid viscosity), plotted as solid symbols to the left side of the plot. Considering $K_I = K_{eff}(v)$ during propagation, the differences between the values computed with eq. (2) (solid symbols) and those estimated from the numerical simulations (empty symbols) are due to the fact that the numerical simulations account for the fluid viscosity, which affects the fluid pressure profile and, consequently, the crack stress intensity factor. The comparison between them confirms that the larger the velocity, the greater the difference between neglecting (solid symbols) or accounting for (empty symbols) the effect of fluid viscosity.

Finally, in Fig. 5a (dashed line), we display the square root function for $K_{eff}(v)$ obtained by combining $\Delta E_f(v) = 6.14v + 0.09$ Pa·m with eq. (8), using a Poisson's ratio $\nu \approx 0.5$, and the 2-D to 3-D conversion function (eq. 4) with $W/L = 0.56$, an averaged value over all the oil-filled injections. In this plot, $K_{eff}(v)$ is extrapolated to the velocity domain of the air-filled fracture propagation experiments, shown as solid symbols to the right side of Fig. 5a. Note that since air is essentially inviscid, $K_I^{(air)}$ values (obtained with eq. 2) should not be significantly affected by neglecting the fluid viscosity. However, it is evident that the relation for $K_{eff}(v)$ still does not match $K_I^{(air)}$ values. This suggests that
<figure><img src="image_5.png" /><figcaption>**Figure 4.** (a)–(c) are the same as Figs 3a–c, but with the numerical simulations (bold lines) accounting for a velocity-dependent fracture energy $\Delta E_f(v) = 6.14 v + 0.09$ Pa·m.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** (a) Log plot of the effective fracture toughness $K_{\text{eff}}$ versus mean velocity $v_{\text{mean}}$ for air- and oil-filled fracture experiments (solid symbols to the left and right side of the plot, respectively) and for oil-filled fracture simulations (empty symbols). The dashed line represents the square root regression of $K_{\text{eff}}$ as function of $v$, obtained by combining the linear regression for $\Delta E_f(v)$ from numerical simulations (Fig. 3d) with eq. (8) for $v = 0.5$, and C being a 2-D–3-D conversion factor obtained from eq. (4). (b) Semi-log plot of the fracture energy $\Delta E_f$ obtained for oil-filled fracture simulations (empty symbols) and air-filled fracture experiments (solid symbols). The vertical bars represent the errors. For each set of $\Delta E_f$ values (for oil- and air-filled cracks), the dashed lines represent the linear regressions as a function of the mean velocity $v_{\text{mean}}$. Symbols stand for different volumes: squares represent 10 mL injections, triangles are for 30 mL, diamonds for 50 mL and circles represent other volumes.</figcaption></figure>

the linear relation for $\Delta E_f(v)$ also cannot be extrapolated over significantly different ranges of propagation velocities, such as those relative to our oil- and air-filled sets of experiments. The linear trend we computed for $\Delta E_f(v)$ improves our numerical simulations (and their fit with the experiments) when applied to a limited velocity range, but extending it over a largely different velocity range does not appear to be appropriated (the fastest oil-filled crack, $\approx 0.1$ mm s⁻¹, propagates at a velocity that is more than one order of magnitude smaller than the slowest air-filled crack, $\approx 2$ mm s⁻¹).

In order to show this, we calculated the effective fracture energy, $\Delta E_f(v)$, for the air-filled crack experiments, first converting $K_f^{3D(air)}$ to 2-D values using eq. (4) (with W/L values relative to each air injection), and then applying eq. (8). We plotted these values alongside those from numerical simulations on Fig. 5b, also displaying the best-fitting linear regressions for both data sets: for the air-filled fracture experiments, we obtained $\Delta E_f^{(air)} = 0.13v + 0.82$ Pa·m (dashed line interpolating solid symbols), while for the oil-filled fracture simulations $\Delta E_f^{(oil)} = 6.14v + 0.09$ Pa·m (dashed line interpolating empty symbols). These linear regressions represent an empirical relationship for the ‘effective fracture energy’ of our gelatin blocks across various velocities and velocity ranges. Notably, at the lower velocities that we tested (the oil-filled fracture experiments and their simulations), the slope of the regression line is significantly steeper than that for the air-filled fracture experiments.

# 4. DISCUSSION

## 4.1. Effective fracture toughness and fracture energy

Our results revealed a significant discrepancy in fracture toughness estimates derived from oil- ($K_c^{(oil)} = 2.2$ Pa·m¹/²) and air-filled crack experiments ($K_c^{(air)} = 24.0$ Pa·m¹/²). These estimates were obtained by extrapolating $K_I = K_{\text{eff}}(v)$ for a vanishing fracture propagation velocity using linear regressions (Fig. 2), applied separately to our data sets of air- and oil-filled cracks. Despite two distinct linear regressions providing a good fit for $K_{\text{eff}}(v)$ values within their respective velocity ranges, the $K_c$ estimates differed by an order of magnitude. To account for the effect of fluid viscosity, we used a numerical model for fluid-filled fracture propagation applied to the oil experiments. The numerical simulations provided a different relationship for $K_{\text{eff}}(v)$ (a square root function), derived by combining a linear regression of $\Delta E_f(v)$ with eq. (8). However, even this square root function for $K_{\text{eff}}(v)$, obtained from oil-filled crack simulations, did not match $K_{\text{eff}}$ values in the velocity range relative to the air-filled crack experiments (cf. Fig. 5a), suggesting that also this relationship for $K_{\text{eff}}(v)$ cannot be extrapolated over several orders of magnitude of velocity.

In fact, one of the main differences between air- and oil-filled fracture propagation experiments—beyond the effects of oil viscosity—is the velocity range of the fracturing process, as air-filled fractures propagate at velocities between 1 and 3 orders of magnitude larger than oil-filled fractures. Additionally, $K_{\text{eff}}$ estimates may depend not
only on fracture propagation velocities but also on fracture lengths, as already suggested by previous studies on hydraulic fracturing (e.g. Papanastasiou 1999), and magmatic dykes (Rubin 1993; Rivalta et al. 2015, and references therein). However, in our data sets, a dependence of $K_{\text{eff}}$ on crack length might be obscured by the more pronounced effect of velocity, given that our crack lengths vary over roughly one order of magnitude (cf. L values in Tables 1 and 2), whereas velocities span almost four (Fig. 5c).

But how relevant are these results for crustal rocks and magmatic intrusions? The large differences that have been observed between fracture toughness values measured from rock samples in the lab and magmatic dykes in the field may be due to a similar effect, which already led Rubin (1993) to suggest that fracture toughness is unlikely a characteristic purely dependent on the rock properties. Laboratory measurements have empirically shown that fracture size may significantly influence fracture toughness. For instance, Ayatollahi & Akbardoost (2014) performed several fracture tests on rock disks with an initial circular crack sized 0.3 times the rock sample diameter. They found that $K_c$ measurements highly depend on the rock sample diameter and, therefore, the initial crack size. Additionally, Rivalta et al. (2015) calculated the ‘critical’ dimensions for the propagation of a magmatic dyke using a field-based estimate for the effective fracture toughness of the crust ($K_{\text{eff}} \approx 100 \text{ MPa} \cdot \text{m}^{1/2}$), which is 2 orders of magnitude larger than values from laboratory measurements, resulting in a critical length and thickness in the order of 2 km and 0.5 m, respectively. Given that these values align with those observed for magmatic dykes, Rivalta et al. (2015) propose a scale-dependency of fracture toughness on dyke dimensions. Rivalta et al. (2015) also showed that the propagation distance travelled by magmatic intrusions displays a dependency on their velocity, with low viscosity magmas propagating longer distances in shorter time ranges (from a few hours to days) under a $K_{\text{eff}}$ in the order of $\approx 1000 \text{ MPa} \cdot \text{m}^{1/2}$.

The point made for the effective fracture toughness also applies to the energy required to generate a new fracture surface, as $K_c$ and $\Delta E_f$ are linked by eq. (7). Hence, the fracture energy should also be a function of the fracture growth velocity rather than a constant parameter, implying that gelatin should not be considered as a linear brittle elastic material over a large range of crack propagation velocities.

Similar arguments may apply to crustal rocks: Zhu et al. (2022) suggested that in addition to the energy dissipation due to viscous flow, a fracture process zone at the crack tip can explain the orders of magnitude difference between laboratory tests on rock samples and field measurements. They found that a fracture process zone can vary from $10^{-2}$ to $10^{-1}$ m in rock samples, which is expected to scale up to $10^1$ m, for magmatic intrusions. In gelatin experiments, a similar effect may be due to a zone undergoing plastic or viscoelastic deformation in the vicinity of the crack tip. Even though our findings for air- and oil-filled cracks cannot be readily scaled up to rock fracture toughness values, we suggest that a similar procedure to the one followed in this study may be applied to magmatic dykes, and help better constrain effective fracture toughness values for magmatic intrusions, if enough observations can be gathered.

We showed that two distinct linear regressions fit our $\Delta E_f$ values at best, within different velocity ranges (e.g. the slower oil-filled crack simulations and the faster air-filled crack injections), displaying an increased $\Delta E_f(v)$ slope at lower velocities (Fig. 5). This suggests that a more general and comprehensive relationship linking $\Delta E_f$ and $v$ should also be able to reproduce such a change in slope. However, to better understand the relationship between the

two observed trends of $\Delta E_f(v)$, new experiments using a fluid with intermediate viscosity are necessary to bridge the gap between the two velocity ranges covered by our sets of experiments. Additional experiments may also highlight the effect of other crack characteristics on $\Delta E_f$, such as crack lengths and overpressure gradients (buoyancy).

## 4.2. Experimental errors and differences between experiments and numerical simulations

We estimated the gelatin fracture energy $\Delta E_f$ at different propagation velocities using a numerical model and fitting the velocity profiles recorded during oil-filled crack propagation experiments with three different volumes: 10, 30 and 50 mL. For the larger volumes (30 and 50 mL), the crack propagation velocity during the initial phase of the experiment may be more affected by the oil injection rate than the 10 mL injections. Indeed, the duration and the velocity profile during the initial crack growth phase depend on the time required to manually inject the viscous oil through the needle. This resulted in a greater dispersion among the velocity profiles relative to the larger volume injections, compared to others (cf. 50 mL velocity profiles, in Fig. 3c).

In order to highlight the impact of this effect on our experimental velocity profiles, we used a thicker needle in experiments 14-O1, 14-O2 and 14-O3, drastically reducing the time needed to inject the oil. By comparing the velocity profiles for experiments 14-O1, 14-O2 and 14-O3 with the others, it is possible to appreciate a larger initial velocity and a faster velocity decay during the injection phase of the experiments performed with a thicker needle, particularly for 30 and 50 mL injections (Fig. S1, Supporting Information). Hence, to ensure that our results are not influenced by the initial velocity decay, we excluded the initial part of the propagation for the fit with the numerical simulations (vertical dotted lines in Fig. 3). Additionally, we noted another difference in the experiments performed with thicker needles: oil-filled cracks developed under higher injection rates, particularly 30 and 50 mL injections, consistently had larger aspect ratios (W/L), mainly due to slightly shorter head lengths L, compared to the cracks developed under slower injection rates (cf. W and L values for exp 14-O1-O2-O3, compared to other experiments with the same volumes). This observation is compatible with the fact that cracks formed under a slower injection rate may have developed in height (due to buoyancy) for a longer time before the injection was finished. Therefore, they accommodated the oil volume along a longer head length with respect to the oil-filled cracks that formed more rapidly under a higher injection rate. However, this difference did not significantly impact velocity values after the initial ‘crack formation phase’, at least within our error estimates (cf. Fig. S1, Supporting Information, and Table 1).

We used eq. (8) to calculate $\Delta E_f$ values for air-filled crack experiments, assuming that the crack propagation occurs with $K_f = K_{\text{eff}}(v)$ (Dahm 2000b), and the viscous dissipation is negligible for air-filled fracture propagation. We also used eq. (8) to compute $K_f$ given the $\Delta E_f$ values obtained from numerical simulations of oil-filled cracks. In both cases, we must consider that values for air-filled crack experiments are in 3-D and numerical simulations for oil-filled cracks are in 2-D. In order to account for the difference introduced by the 2-D approximation of the numerical model and compare numerical simulations and experiments, we used eq. (4) with measured values for W/L.

In particular, we used the averaged W/L values for 10, 30 and 50 mL oil injections when computing 3-D estimates of $K_{\text{eff}}$ based
<table><caption>**Table 3.** Average parameter values used to constrain the numerical simulations, their standard deviation (assuming that they follow a Gaussian distribution), and the ratio between their average and standard deviation.</caption><thead><tr><th>Vol</th><th>ID</th><th>ρ</th><th>σ<sub>ρ</sub></th><th>σ<sub>ρ</sub>/ρ [%]</th><th>G</th><th>σ<sub>G</sub></th><th>σ<sub>G</sub>/G [%]</th><th>O</th><th>σ<sub>O</sub></th><th>σ<sub>O</sub>/O [%]</th></tr></thead><tbody><tr><td>10</td><td>10-AV</td><td>1110</td><td>10.6</td><td>0.96</td><td>135</td><td>7.3</td><td>5.4</td><td>0.45</td><td>0.06</td><td>13.3</td></tr><tr><td>30</td><td>30-AV</td><td>1107</td><td>1.5</td><td>0.13</td><td>133</td><td>8.0</td><td>6.0</td><td>0.78</td><td>0.07</td><td>9.0</td></tr><tr><td>50</td><td>50-AV</td><td>1106</td><td>3.0</td><td>0.27</td><td>134</td><td>8.3</td><td>6.2</td><td>0.98</td><td>0.03</td><td>3.1</td></tr></tbody></table>

on $\Delta E_f$ values from 2-D numerical simulations. Similarly, to provide 2-D estimates for $\Delta E_f$ based on 3-D values of $K_{eff}$ from air-filled crack experiments, we used the $W/L$ ratios from each of those experiments. Finally, for plotting a 3-D estimate for the effective fracture toughness $K_{eff}(v) = \sqrt{4G} \cdot (6.14v + 0.09)$ derived from the 2-D numerical simulations, we used the average $W/L$ ratio over all oil-filled cracks. Such average value does not capture the trends exhibited by $W/L$ ratios as a function of $L$, shown in the Supporting Information, Fig. S2, where we plot $W(L)$ for both sets of experiments. These different $W(L)$ trends for oil- and air-filled cracks are not trivial to interpret and may suggest additional dependencies on the crack geometry that have not been fully considered in our analysis, such as effects related to injection rates. However, within the range of variability of $W/L$ values, our results remain robust, as evidenced by the small differences between the dashed line and symbols in Fig. 5a.

As we mentioned in Section 2.1.2, errors associated with length measurements may vary according to the sharpness and possible distortion effects of the recorded images. Particularly, since the camera used to capture crack opening and width was not fixed to the experimental setup, the reference scale sometimes suffered from distortion effects that increased measurement errors. This mainly affected the errors in $K_I$, which are strongly influenced by inaccuracies in crack length and width measurements (eq. 2). In general, the relative errors associated with parameter estimates can range from a few per cent to approximately 50%, or—in exceptional cases—even higher, particularly for the lowest values of propagation velocities (Table 1). These relatively large errors primarily stem from manually picking the crack tip positions and the sharpness of the images used for this purpose.

Smittarello et al. (2021) compared dynamic Young’s modulus estimates obtained with shear wave velocity measurements in gelatin with static estimates obtained with deformation measurements induced by a load placed at the gelatin’s surface. Their results showed that Young’s modulus estimates with these two methods may sometimes display large differences, but they could not highlight any systematic trend. In their study, the measurements were performed on gelatin with different concentrations, making it difficult to understand which of the two methods provided more stable results. On the other hand, our experiments are performed using always the same gelatin concentration, and the fact that our shear modulus estimates are consistent (within errors) indicates that the technique we used provided robust estimates. Smittarello et al. (2021) also noted that the shear wave method reflects the gelatin condition within the block. In contrast, estimates obtained using the surface load method only provide information about the gelatin’s shallow portion, failing to offer insights into the gelatin’s condition at greater depths. Additionally, surface load estimates may be significantly influenced by the ratio of the load surface to the distance from the walls, and they exhibit stiffness gradients depending on where the load is placed. For these reasons, we consider the shear wave method preferable for our modelling purpose. Eventually, we note that the shear modulus of our gelatin blocks (≈ 100 Pa) would correspond to

approximately 100 MPa in nature (cf. Kavanagh et al. 2013), which may be 1 to 2 orders of magnitude lower than crustal rocks. However, increasing the gelatin rigidity would also increase its strength, making oil-filled crack propagation more difficult, requiring larger crack lengths that would not be manageable with our laboratory setup.

In order to constrain the numerical simulations, we computed average values for several parameters relative to the oil-filled fracture experiments for 10, 30 and 50 mL injections (10-AV, 30-AV and 50-AV rows in Table 1). The uncertainties associated with these average parameters were determined by taking the greater of the two values: (1) the standard deviation of the parameter measurements for each set of experiments with the same injected volume, and (2) the average measurement error for each set of experiments with the same injected volume. This approach ensures that the reported uncertainty reflects the actual variability of each parameter across different experiments, as indicated by the standard deviation unless this variability is smaller than the average instrumental error associated with those measurements. This conservative approach was chosen because a low standard deviation may sometimes be merely due to the relatively small number of experiments performed with the same injected volume.

Eventually, in order to provide uncertainties on the $\Delta E_f$ values obtained from the numerical simulations, we performed additional simulations constrained by the values of the maximum fracture opening plus and minus their associated uncertainties ($O$ and $\sigma_O$ in 10-AV, 30-AV and 50-AV rows of Table 1). Among the input parameters for the numerical simulations—gelatin and oil densities, gelatin rigidity and maximum fracture opening ($O$)—we chose to vary $O$, as it displayed the larger variability among experiments with the same injected volumes (Table 3).

## 5. CONCLUSION

We used two sets of experiments with two different fluid viscosities to estimate the fracture toughness of gelatin blocks with the same elastic properties. We first provided these fracture toughness estimates using static crack formulae that neglect the effect of fluid viscosity (e.g. Weertman 1971). Through this approach, we found that fracture toughness estimates obtained from air- and oil-filled crack experiments differ by more than 1 order of magnitude. The equation used to derive these estimates does not account for the (viscosity-dependent) dynamic fluid-pressure profile, making it less suitable for high-viscosity oil-filled cracks. To address this limitation, we used the numerical model from Furst et al. (2023), which accounts for the effect of fluid viscosity and allowed us to provide new estimates of the crack stress intensity factor based on the simulations of three different volumes of oil-filled crack experiments. Our results show that when the fluid viscosity is considered, the faster the crack, the larger the difference between stress intensity factor estimates that neglect or account for the fluid viscosity (solid and empty symbols to the left side of Fig. 5a). Furthermore, our
analysis demonstrates that the estimates of $\Delta E_f(v)$ derived from the numerical simulations of oil-filled experiments and those from air-filled experiments cannot be adequately described by a single linear regression across a velocity range spanning several orders of magnitude (Fig. 5b). This may be partly attributed to the fracture energy and, therefore, the effective fracture toughness of gelatin, which may not only depend on the propagation velocity but also on other crack characteristics such as its dimension and/or internal pressure gradient (buoyancy). These findings may be very relevant for magmatic dykes, with similar arguments for the effective fracture toughness applying to crustal rocks, as it has been previously proposed and discussed by several authors (e.g. Heimpel & Olson 1994; Rubin 1998; Rivalta et al. 2015; Zhu et al. 2022). We conclude that introducing an effective, velocity-dependent fracture toughness in theoretical frameworks is necessary when aiming to describe fluid-filled fracture propagation velocities across various regimes. The main challenge lies in establishing empirical relationships for $K_{\text{eff}}(v)$ at the scale of magmatic intrusions, similar to our approach in analogue experiments. Achieving this would require extensive direct observations of the propagation velocity and characteristics of magmatic dykes within crustal rocks. Such comprehensive data collection would necessarily involve monitoring observations from several active volcanic systems.

## FUNDING

This research was funded by the INGV project 'BEST-CaFE' (Pianeta Dinamico). SF and VP were financially supported by the ANR-DFG NLE 2018 MagmaPropagator project (ANR-18-CE92-0037). FM was supported by the VOLC-LABS project (BRL2019.19, funded by INGV).

## ACKNOWLEDGMENTS

Thanks to the collaboration with the German Research Centre for Geosciences, the experiments were performed at the 'Volcano Lab', GFZ—Section 2.1, Potsdam. We especially thank Lorenzo Mantiloni for helping carry out the lab models, and we thank the lab technicians who contributed to building up the experimental setup: Ralf Bauz and Stefan Mikulla. We also thank the editor Shiqing Xu, Andreas Möri and Thierry Menand for their comments and suggestions, which significantly contributed to improving this article.

## SUPPORTING INFORMATION

Supplementary data are available at *GJIRAS* online.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper.

## DATA AVAILABILITY

Experimental data and code used for simulations are available on the repository: https://zenodo.org/record/8271248.

## REFERENCES

* Ayatollahi, M.R. & Akbardoost, J., 2014. Size and geometry effects on rock fracture toughness: Mode I fracture. *Rock Mech. Rock Eng.*, **47**, 677–687.

* Bhat, H.S., Rosakis, A.J. & Sammis, C.G., 2012. A micromechanics Based-Constitutive model for brittle failure at high strain rates. *J. appl. Mech.*, **79**, 031016, doi:10.1115/1.4005897

* Brizzi, S., Funiciello, F., Corbi, F., Giuseppe, E. Di & Mojoli, G., 2016. Salt matters: how salt affects the rheological and physical properties of gelatine for analogue modelling. *Tectonophysics*, **679**, 88–101.

* Brown, D., 2012. Tracker video analysis and modeling tool for physics education. Retrieved September 20, 2017, from http://physlets.org/track er/

* Chen, L., Zhang, G., Zou, Z., Guo, Y. & Zheng, X., 2021. The ef fect of fracture growth rate on fracture process zone development in quasi-brittle rock. *Eng. Fract. Mech.*, **258**, 108086, Pergamon, doi:10.1016/J.ENGFRACMECH.2021.108086

* Corvec, N. Le, Menand, T. & Lindsay, J., 2013. Interaction of ascending magma with pre-existing crustal fractures in monogenetic basaltic volcanism: an experimental approach. *J. geophys. Res. Solid Earth*, **118**, 968–984.

* Dahm, T., 2000a. On the shape and velocity of fluid-filled fractures in the Earth. *Geophys. J. Int.*, **142**, 181–192.

* Dahm, T., 2000b. Numerical simulations of the propagation path and the arrest of fluid-filled fractures in the Earth. *Geophys. J. Int.*, **141**, 623–638.

* Davis, T., Rivalta, E., Smittarello, D. & Katz, R.F., 2023. Ascent rates of 3-D fractures driven by a finite batch of buoyant fluid. *J. Fluid Mech.*, **954**, A12, doi:10.1017/JFM.2022.986

* Delaney, P.T. & Pollard, D.D., 1981. Deformation of Host Rocks and Flow of Magma during Growth of Minette Dikes and Breccia-bearing Intrusions near Ship Rock, New Mexico, USGPO. **1202**. https://doi.org/10.3133/pp 1202.

* Delaney, P.T., Pollard, D.D., Ziony, J.I. & McKee, E.H., 1986. Field relations between dikes and joints: emplacement processes and paleostress analysis. *J. geophys. Res. Solid Earth*, **91**, 4920–4938.

* Furst, S., Maccaferri, F. & Pinel, V., 2023. Modeling the shape and ve locity of magmatic intrusions, a new numerical approach. *J. geophys. Res. Solid Earth*, **128**, e2022JB025697, John Wiley & Sons, Ltd, doi: 10.1029/2022JB025697

* Furst, Séverine, Pinel, V. & Maccaferri, F., 2024. Dynamics of magmatic intrusion: what can we learn from the comparison of analog and numerical models? *Volcanica*, **7**, 67–87.

* Garagash, D. & Germanovich, L., 2014. Gravity driven hydraulic fracture with finite breadth, In *Proc. Soc. Eng. Sci. 51st Annual Technical Meet ing*, eds. Bajaj, A., Zavattieri, P., Koslowski, M. & Siegmund, T., West Lafayette: Purdue University Libraries Scholarly Publishing Services. https://docs.lib.purdue.edu/ses2014/honors/rudnicki/13

* Garagash, D. & Germanovich, L., 2022. Notes on propagation of 3D buoyant fluid-driven cracks. arXiv preprint arXiv:2208.14629

* Heimpel, M. & Olson, P., 1994. Buoyancy-driven fracture and magma trans port through the lithosphere: models and experiments. *Int. Geophys.*, **57**, Academic Press Inc, doi: 10.1016/S0074-6142(09)60098-X

* Kavanagh, J.L., Menand, T. & Daniels, K.A., 2013. Gelatine as a crustal ana logue: determining elastic properties for modelling magmatic intrusions. *Tectonophysics*, **582**, 101–111, Elsevier B.V.

* Kavanagh, J.L., Engwell, S.L. & Martin, S.A., 2018. A review of laboratory and numerical modelling in volcanology. *Solid Earth*, **9**, 531–571.

* Lengliné, O., Duputel, Z. & Okubo, P.G., 2021. Tracking dike propagation leading to the 2018 Kīlauea eruption. *Earth planet. Sci. Lett.*, **553**, 116653, Elsevier, doi:10.1016/J.EPSL.2020.116653

* Lister, J.R., 1990. Buoyancy-driven fluid fracture: the effects of material toughness and of low-viscosity precursors. *J. Fluid Mech.*, **210**, 263–280.

* Liu, D. & Lecampion, B., 2022. Laboratory investigation of hydraulic frac ture growth in Zimbabwe Gabbro. *J. geophys. Res. Solid Earth*, **127**, 127, doi:10.1029/2022JB025678

* Liu, D. & Lu, G., 2023. Effects of velocity-dependent apparent toughness on the pre- and post-shut-in growth of a hydraulic fracture. *Comput. Geotech.*, **155**, 105195, Elsevier Ltd, doi:10.1016/j.compgeo.2022.105195

* Maccaferri, F., Bonafede, M. & Rivalta, E., 2010. A numerical model of dyke propagation in layered elastic media. *Geophys. J. Int.*, **180**, 1107–1123.
Möri, A. & Lecampion, B., 2022. Three-dimensional buoyant hydraulic fractures: constant release from a point source. *J. Fluid Mech*, **950**, 12, doi:10.1017/jfm.2022.800

Möri, A. & Lecampion, B., 2023. Three-dimensional buoyant hydraulic fractures: finite-volume release. *J. Fluid Mech.*, **972**, 20, doi:10.1017/jfm.2023.711

Olson, J.E. & Schultz, R.A., 2011. Comment on ``A note on the scaling relations for opening mode fractures in rock'' by C.H. Scholz, *J. Struct. Geol.*, **33**, 1523-1524.

Pansino, S. & Taisne, B., 2019. How magmatic storage regions attract and repel propagating dikes. *J. geophys. Res. Solid Earth*, **124**, 274-290.

Papanastasiou, P., 1999. The effective fracture toughness in hydraulic fracturing. *Int. J. Fract.*, **96**, 127-147.

Pinel, V., Furst, S., Maccaferri, F. & Smittarello, D., 2022. Buoyancy versus local stress field control on the velocity of magma propagation: insight from analog and numerical modelling. *Front. Earth Sci. (Lausanne)*, **10**, Frontiers Media SA, doi:10.3389/feart.2022.838318

Prejean, S., Stork, A., Ellsworth, W., Hill, D. & Julian, B., 2003. High precision earthquake locations reveal seismogenic structure beneath Mammoth Mountain, *Geophys. Res. Lett.*, **30**, Wiley-Blackwell, doi:10.1029/2003GL018334

Rice, J.R., 1968. A path independent integral and the approximate analysis of strain concentration by notches and cracks. *J. appl. Mech.*, **35**, 379-386.

Rivalta, E., Taisne, B., Bunger, A.P. & Katz, R.F., 2015. A review of mechanical models of dike propagation: schools of thought, results and future directions. *Tectonophysics.*, **638**, 1-42.

Rivalta, E. & Dahm, T., 2006. Acceleration of buoyancy-driven fractures and magmatic dikes beneath the free surface. *Geophys. J. Int.*, **166**, 1424-1439.

Roper, S.M. & Lister, J.R., 2007. Buoyancy-driven crack propagation: the limit of large fracture toughness. *J. Fluid Mech.*, **580**, 359-380.

Rubin, A.M., 1993. Tensile fracture of rock at high confining pressure: implications for dike propagation. *J. geophys. Res. Solid Earth*, **98**, 15919-15935.

Rubin, A.M., 1995. Propagation of magma-filled cracks. *Annu. Rev. Earth Planet. Sci.*, **23**, 287-336. Retrieved from https://www.annualreviews.org/doi/pdf/10.1146/annurev.ea.23.050195.001443

Rubin, A.M., 1998. Dike ascent in partially molten rock. *J. geophys. Res.—Solid Earth*, **103**, 20901-20919.

Secor, D.T. & Pollard, D.D., 1975. On the stability of open hydraulic fractures in the Earth's crust. *Geophys. Res. Lett.*, **2**, 510-513.

Sigmundsson, F., et al., 2015. Segmented lateral dyke growth in a rifting event at Bárðarbunga volcanic system, Iceland. *Nature*, **517**, 191-195.

Sigmundsson, F., et al., 2024. Fracturing and tectonic stress drives ultrarapid magma flow into dikes, *Science* (1979), **383**, 1228-1235.

Smittarello, D., Pinel, V., Maccaferri, F., Furst, S., Rivalta, E. & Cayol, V., 2021. Characterizing the physical properties of gelatin, a classic analog for the brittle elastic crust, insight from numerical modeling. *Tectonophysics*, **812**, 228901, doi:10.1016/j.tecto.2021.228901

Spence, D.A., Sharp, P.W. & Turcotte, D.L., 1987. Buoyancy-driven crack propagation: a mechanism for magma migration. *J. Fluid Mech.*, **174**, 135-153.

Spence, D.A. & Turcotte, D.L., 1990. Buoyancy-driven magma fracture: A mechanism for ascent through the lithosphere and the emplacement of diamonds. *J. geophys. Res. Solid Earth*, **95**, 5133-5139.

Takada, A., 1990. Experimental study on propagation of liquid-filled crack in gelatin: shape and velocity in hydrostatic stress condition. *J. geophys. Res.*, **95**, 8471-8481.

Watanabe, T., Masuyama, T., Nagaoka, K. & Tahara, T., 2002. Analog experiments on magma-filled cracks: competition between external stresses and internal pressure. *Earth Planets Space*, **54**, e1247-e1261.

Weertman, J., 1971. Theory of water-filled crevasses in glaciers applied to vertical magma transport beneath oceanic ridges, *J. geophys. Res.*, **76**, 1171-1183.

Zhu, D., Han, G., Zou, H., Cui, M., Liang, C. & Yao, F., 2022. A review of the hydraulic fracturing in ductile reservoirs: theory, simulation, and experiment, *Processes*, **10**, 2022, doi:10.3390/PR10102022
# Supporting Information for "Investigating the velocity of magmatic intrusions and its relation with rock fracture toughness: insights from laboratory experiments and numerical models"

A. Gaete<fnref n="1" />,<fnref n="2" />,<fnref n="3" />, F. Maccaferri <fnref n="1" />, S. Furst<fnref n="4" />, V. Pinel<fnref n="5" />

## Contents of this file

* Table S1: 3D volume estimates for simulations

* Figure S1: Experimental velocity profiles

* Figure S2: Experimental W/L ratio for air- and oil-filled cracks

## Introduction

This supporting information provides extra material with the table of 3D volume estimates for the simulated cracks, the plots -separated by volumes- of the velocity profiles tracked from the experiments and used for simulations, and the linear regression results that demonstrate the relationship between the width (W) and length (L) of the crack.

<footnote n="1">¹ INGV Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Vesuviano, Naples, Italy</footnote>

<footnote n="2">² GFZ German Research Centre for Geosciences, Potsdam, Germany</footnote>

<footnote n="3">³ OVDAS Observatorio Volcanológico de Los Andes del Sur, SERNAGEOMIN, Temuco, Chile</footnote>

<footnote n="4">⁴ GEOMAR Helmholtz-Centre for Ocean Research, Kiel, Germany</footnote>

<footnote n="5">⁵ University Grenoble Alpes, University Savoie Mont Blanc, CNRS, IRD, University Gustave Eiffel, Grenoble, France</footnote>
<table><caption>**Table 1.** 3D volume estimates for the simulations, obtained by multiply $A_0$ times the average width for each volume experiment.</caption><thead><tr><th>Injected volumes (mL)</th><th>A<sub>0</sub> (cm<sup>2</sup>)</th><th>Mean width (cm)</th><th>~ 3D volume estimate (mL)</th></tr></thead><tbody><tr><td>10</td><td>2.4</td><td>3.4</td><td>8</td></tr><tr><td>30</td><td>4.9</td><td>5.6</td><td>27</td></tr><tr><td>50</td><td>6.7</td><td>6.6</td><td>44</td></tr></tbody></table>
<figure><img src="image_7.png" /><figcaption>**Figure 1.** Experimental velocity profiles. Velocity measured along the propagation of oil-filled cracks for a) 10 mL, b) 30 mL and c) 50 mL injections. The black curves represent experiments performed with a thicker needle: 14-O1, 14-O2 and 14-O3 for 10 mL, 30 mL and 50 mL, respectively (Table 1).</figcaption></figure>
<figure><img src="image_8.png" /></figure>

**Figure 2.** Experimental W/L ratio for a) air- and b) oil-filled cracks. The graphs depict the crack width (W) in proportion to its length (L) for a comprehensive range of injections (as detailed in Table 2 and 1 respectively). The dashed line represents the optimal linear regression, with a formula of $f(x) = 1.31x + 1.59$ for air-filled cracks. Meanwhile, the solid line demonstrates a linear regression with a constraint on the origin of the coordinates. This suggests that as the volume decreases towards zero, both W and L must approach zero, as observed in Dahm (2000a). Although the dashed line offers the best fit for the experimental data, it is not physically plausible because it does not intersect the (W,L) plane at the origin (0,0). Note the difference in scale between the two subfigures.