Evaluation of the submarine debris-flow hazard risks to planned subsea pipeline systems: a case study in the Qiongdongnan Basin, South China Sea

1. Introduction

China is a major consumer of oil and gas resources. In 2021, China’s dependence on foreign oil and natural gas exceeded 70% and 45%, respectively, and this dependence has continued to rise (Liu et al., 2022). Fortunately, China is rich in oil and gas resources in the deepwater area of the South China Sea. As such, there is a strategic need to guarantee national energy security by accelerating the development of deepwater oil and gas fields in the South China Sea. The successful development and production of the Liwan 3-1 deepwater gas field (Jin et al., 2018), Liuhua 16-2 deepwater oil field (Wang et al., 2020), and Lingshui 17-2 deepwater gas field (Zhu et al., 2018) in China have kicked off the development of deepwater oil and gas fields in the South China Sea, and more new developments are planned in the future.

One recent development is currently underway in the South China Sea. This new project is located approximately 132 km southwest of Sanya City, Hainan Island, which stretches across the continental shelf edge and parts of the continental slope down to an approximate water depth of 890 m to 970 m. This area is part of the Qiongdongnan Basin, a major petroliferous basin in the northern South China Sea, extending over an area of approximately 45 000 km². The Qiongdongnan Basin is bound by Hainan Island in the north, the Zhujiang River Mouth Basin to the east, the Xisha Uplift to the south, and the Yinggehai Basin and Red River Fault System to the west (Su et al., 2016). There is significant evidence that this area has been subjected to repeated mass transport events, typically in the form of debris gravity flows, during periods of its geological history. The mass transport deposits (MTDs) in the representative form of submarine debris flows occurred because of high sedimentation rates since the late Miocene and earthquake activities (Xie et al., 2008; Wu et al., 2009; Li et al., 2015; Wang et al., 2019; Cheng et al., 2021).

Submarine debris flows are not only an important method of massive sediment transport from the continental shelf to the deepwater basin, but also one of the severe marine geological hazards (Li et al., 2012; Xiu et al., 2015; Jia et al., 2016; Wang et al., 2016, 2019; Kim et al., 2019; Nugraha et al., 2022). Upon initiation of submarine landslides, the material strength of seafloor sediments with high sensitivity significantly decreases; meanwhile the bulk of released sediments might quickly mix with the ambient seawater during downslope sliding processes (Elverhøi et al., 2005, 2010; Ilstad et al., 2004a, b; Zakeri et al., 2008; Li et al., 2012; Haza et al., 2013; Xiu et al., 2015; Liu et al., 2018; Guo et al. 2022a; Du et al., 2022). As a result, the failed materials tend to rapidly transform into the submarine debris flows. Driven by gravity, they can be in motion with very high velocities for less than an hour to several days, and their runout distances might vary significantly from less than several tens of meters to hundreds of kilometers on gentle seabed (Imran et al., 2001; Drago, 2002; Niedoroda et al., 2003; Gauer et al., 2005; De Blasio et al., 2005; Elverhøi et al., 2005; Das, 2012; Xiu et al., 2015; Kim et al., 2019). Due to the powerful momentum and vast sphere of influence, submarine debris flows have become one of the most notorious marine geohazards that can potentially devastate various seafloor infrastructures, especially long-spanned oil and gas pipelines installed across their runout path, which has become an increasingly operative issue in the offshore pipeline engineering industry (Niedoroda et al., 2003; Bruschi et al., 2006; Zakeri et al., 2008; Das, 2012; Haza et al., 2013; White et al., 2016; Malgesini et al., 2018; Nian et al., 2018; Chaytor et al., 2020; Xiu et al., 2021; Guo et al., 2022a).

For a new deepwater oil and gas development in the Qiongdongnan Basin, submarine pipelines are one of the key components of the whole development and production system, and they are indispensable for the safe transportation of large quantities of oil and gas from seabed wells to land terminals (Bruschi et al., 2006; Zakeri et al., 2008; Haza et al., 2013; White et al., 2016; Malgesini et al., 2018; Chaytor et al., 2020; Xiu et al., 2021). A preliminary pipeline alignment was proposed across the targeted development region. However, the pipeline route was planned to traverse this specific development region where geo-hazard risks arising from future potential submarine debris flows might occur. It is necessary to take this geological risk into account when planning a pipeline route because such hazardous mass transport events have been observed to occur historically during geological time scales (Xie et al., 2008; Wu et al., 2009; Li et al., 2015; Wang et al., 2019; Cheng et al., 2021). Once such future potential debris-flow events occur and run across the pipeline alignment, they would likely cause incalculable losses in economic and environmental terms. As a result, there is a pressing need to evaluate the debris-flow hazard risks to planned pipeline systems in this development region and provide scientific support to design a more resilient pipeline alignment in consideration of future debris-flow hazard risks.

To this end, a numerical code (Spinewine et al., 2011) developed in-house by Fugro was used to perform a series of simulations to support the quantitative evaluation of debris-flow hazard risks to proposed submarine pipeline systems in the development site within the Qiongdongnan Basin, South China Sea. The work started with the back-analysis of a historical debris-flow event and the calibration of a failure material runout model. Then the numerical code was used to predict influence zones of the debris flows originating from potential future slope failures. The spatiotemporal variations of debris-flow thickness and velocity as they propagated down across the natural seafloor topography were also obtained. The modeling results of this case study were adopted during the quantitative evaluation of the hazard risks imposed to submarine pipeline systems located within the influence zone of runout debris materials, and to determine whether the planned pipeline alignment should be adjusted to minimize the hazard risks of potential pipeline failures in such events.

The rest of this paper is organized as follows. In Section 2, the data and methodology used in this quantitative evaluation are briefly described. In Section 3, the modeling results and analyses are presented. In Section 4, the limitations of the modeling results are discussed, and two alternatives for the optimization of the proposed pipeline routes are recommended. In Section 5, the conclusions are drawn.

2. Data and methodology

2.1 Basic datasets

2.1.1   Proposed infrastructure layout

Several pipelines have been proposed at this new development site. These pipelines connect three Manifolds A, B, and C, and they cross the upper continental slope at a less steep section to shallow shelf waters. The planned alignment of the three manifolds and associated pipeline systems is shown (Fig. 1). Note, the fictitious layout of pipelines and manifolds is shown instead of the specific locations of these seafloor infrastructures.

Figure  1.  Planned pipeline systems and mapped mass transport deposits (MTD) within the development site.

DownLoad: Full-Size Img PowerPoint

2.1.2   Bathymetric data

The bathymetric data consisted of the measured depth from the sea surface to the seabed. Four datasets were acquired by a multi-beam echo sounder across the development region, and three datasets were derived from the seafloor picks of the available three-dimensional (3D) seismic data. Each dataset had different spatial extents to cover the targeted areas, and lateral resolutions ranging from 1 m to 15 m bin sizes. These data were integrated into a single bathymetry surface while maintaining the resolution as high as possible in each region. The composite bathymetric surface was used in this modeling study (Fig. 1).

2.1.3   3D seismic data

Historical MTDs can be mapped from 3D seismic data. In this study, the three 3D seismic data sets were available across the development area. However, no relic failure was identified at the seabed surface, as drape sediments with tens of meters of thickness covered the entire study area. Only one of the recent MTDs was mapped from the 3D seismic data. The MTD as well as the example inline (Line 7343) and crossline (Crossline 6628) from 3D seismic data are presented (Fig. 2a). The mapped MTD lay several tens of meters below the seabed, and it was covered by several tens of meters of sediment drape. It was estimated that the age of mapped MTD was at least 200 000 years given its depth below the present-day seafloor. Given that no age dating had been performed on the sediments, this age was estimated by applying an average sedimentation rate of 455 m/Ma since the deposition of the MTD. It was uncertain that whether the interpreted MTD was composed of multiple but smaller individual MTDs, that is, successive events over a short period of the geological time, or represented one large deposit. The uncertainty was due to limited vertical resolution of the 3D seismic data (Fig. 2b). The MTD likely continued significantly beyond interpreted extent with reduced thickness. However, this smaller thickness below 8 m could not be resolved in the available 3D seismic data. There was evidence in some of the 3D seismic lines for boundaries that can divide the mapped MTD into two or more individual MTDs. However, these boundaries could not be mapped throughout the respective area with enough confidence because they were not consistent across all of the 3D seismic data. Given the uncertainties associated with the mapped deposits, it was assumed that the mapped MTD was a stacked deposit composed of many smaller deposits with thicknesses on the order of 0.5–2 m. It implied much lower material strengths than if the deposit were a single, thick deposit.

Figure  2.  Mapped mass transport deposits (MTD): two example profiles of Line 7343 and Crossline 6628 (a), and uncertainties associated with the mapped deposits (b).

DownLoad: Full-Size Img PowerPoint

2.2 Evaluation methodology

2.2.1   Model descriptions

The Fugro slope stability analyses have shown that potential slope failures were the most likely to occur either during or following a seismic event, based on the understanding of geotechnical conditions in the development site. The submarine debris flows might occur under certain conditions following failures of the slopes, while the shear strains that developed in shallow soils entered the strain-softening range and were not arrested by any immediate reduction in the slope gradient. To quantify such hazard risks to the proposed seafloor infrastructures across the development site from mass movement events developing subsequent to slope failures, the runout behaviors of several possible debris flows originating from the unstable locations were simulated. The runout processes of debris flows are usually modeled using finite difference flow models based on the depth-integrated equations of mass and momentum conservation, which are solved based on either an Eulerian or Lagrangian framework. The representative examples of such approaches include the BING program (Imran et al., 2001), which models the one-dimensional (1D) spreading of debris flow down a slope, and the DM-2D (Niedoroda et al., 2006) and MassFlow2D (Qian et al., 2020) models, which are developed to simulate the two-dimensional (2D) runout of debris flows. In this study, numerical simulations were performed using an in-house code developed by Fugro (Spinewine et al., 2011). The conservation equations of momentum and mass are:

where h is the thickness of the debris flow (m), t is the time (s), $ u $ and $ v $ are the depth-integrated velocities in the x and y directions (m/s), $ \rho $ is the density of the debris flow (kg/m³), $ {\rho _{\text{a}}} $ is the density of the ambient water (kg/m³), $ g $ is the gravitational acceleration (m/s²), $ {z_{\text{b}}} $ is the bed elevation (m), and ${\tau _{{{{\rm{b}}}}{\rm{x}}}}$ and ${\tau _{{{{\rm{b}}}}{\rm{y}}}}$ are the bed shear stresses resisting motion in the x and y directions (Pa).

Submarine debris flows can be commonly assumed to be very dense non-Newtonian fluids in laminar regimes (Imran et al., 2001; Zakeri et al., 2008; Xiu et al., 2015; Guo et al., 2022b). There are many rheological models, such as the Power-law, Bingham, and Herschel–Bulkley models, to describe the rheological behaviors of a non-Newtonian fluid. If shear stress is applied on a Power-law fluid, the associated shear rate appears immediately. However, a yield stress in a Herschel–Bulkley or Bingham fluid must be overcome before the occurrence of shear rate and initiation of fluid motion. Since most experiments show that the debris flow has a clear yield stress, the Herschel–Bulkley and Bingham models are usually preferred to describe the rheological property of debris flows (Imran et al., 2001; Zakeri et al., 2008; Haza et al., 2013; Nian et al., 2018; Qian et al., 2020; Du et al., 2022). Note, the Bingham model is a reduced formulation of the Herschel–Bulkley model with the power index set at unity. The rheological behavior of a submarine debris flow is closely matched with the Herschel–Bulkley model in the full range of shear rate, whereas it can only be fitted based on the Bingham model at relatively high shear rates. Thereby, the Herschel–Bulkley model is selected to describe the rheology of debris flow in this study. The governing equation of Herschel–Bulkley model is expressed as follows (Herschel and Bulkley, 1926):

where $ \tau $ is the shear stress (Pa), $ \dot \gamma $ is the shear strain rate (s⁻¹), $ {\tau _{\text{y}}} $ is the yield stress (Pa), K is the consistency index (Pa·sⁿ), and n is the power index.

The numerical model relied on an unstructured triangular mesh for the representation of the bathymetry. The mesh was constructed using the GMSH, an open-source mesh generator. This mesh generator can provide the capability for local mesh refinements in areas of anticipated runout and along steep slopes and areas with channels to achieve best accuracy in the area of interest. The governing equations were solved using a finite volume approach and a Riemann solver based on the Harten–Lax–van Leer scheme. A lateralized treatment of nonconservative products was used to achieve robustness and accuracy in areas of steep gradients (Spinewine et al., 2011). The solver can be run at first-order and second-order accuracy based on the well-known monotonic upstream-centered scheme for conservation laws gradient interpolation with a minmod limiter and predictor-corrector time integration. However, the in-house developed solver was found to be more stable numerically at first-order accuracy. In addition, the model setting with second-order accuracy was found to be much more time-consuming. As a result, based on the Fugro engineering experiences and recommendations, simulations in this study were performed at first-order accuracy to limit computer processing time constraints and to prevent any potential numerical instabilities across the wetting and drying fronts.

2.2.2   Phases of runout analyses

The in-house numerical code developed by Fugro has been validated thoroughly based on the (1) static conditions with two fluids of identical densities at rest above each other, (2) dam-break flows involving collapse of a body of light fluid over a uniform layer of a denser fluid, and (3) the Liska and Wendroff’s ill-posed test examples (Liska and Wendroff, 1997) involving the two-layer flows over a topographic bump (Spinewine et al., 2011). It has also been applied to perform debris flow runout analyses for the purpose of debris-flow hazard evaluations in many deepwater oil and gas development projects (George Zhang, personal communication, September 20, 2022). In the evaluation activities of debris-flow hazard risks, many back-analyses were performed by Fugro engineers at an attempt to justify the selections of rheological parameters of debris flows for specific sites under consideration. However, this model could not be calibrated against any true debris flow events, as the unknown seafloor conditions prior to past events have prevented any true validation. Despite the shortcomings, the same modeling strategy will still be used in this case study to evaluate the debris-flow hazard risks to planned submarine pipeline systems in the Qiongdongnan Basin, South China Sea. In general, the runout analyses of submarine debris flows detailed in this case study were separated into two phases, which are summarized as follows.

Phase 1 was the back analysis of a historical event. The relict failure visible in geophysics acquired across the targeted development site was modeled numerically, which was based on a set of rheological parameters describing the soil’s resistance to flow. In lieu of site-specific rheology tests, these rheological parameters were selected based on engineering experiences in similar soils in this region. The back-analysis modeling was used to examine the applicability of rheological parameters before using them to forecast runout behaviors of potential future debris flows at the development site.

Phase 2 was the forward prediction analyses of future potential event. This phase included the numerical modeling of debris flows originating from slope stability profiles considered to be most critical to the planned Manifold A, identified through the slope stability analyses. Based on the outcomes of modeling efforts, a similar runout analysis was performed to consider another slope stability profile that was qualitatively considered to be the most critical to seafloor infrastructures surrounding the planned Manifold B. These modeling efforts relied on the rheological parameters settled upon in Phase 1 analysis and engineering experiences in similar soils within this region.

2.3 Modeling specifics

2.3.1   Rheological parameters of debris flows

In the modeling exercise, appropriate rheological parameters had to be selected. As no rheology tests had been conducted on the soil materials recovered from the development site, the rheological parameters were selected based on engineering experiences in similar soils in the region. To develop credible bounds for the runout behaviors of debris flows, three different variations on the material rheology were considered in this study. The Best Estimate (BE) rheological parameters corresponded to the best guess of the rheology of the materials at the development site and were validated in the Phase 1 back analysis. The Low Estimate (LE) and High Estimate (HE) rheological parameters used in the Phase 2 forward analyses were based on engineering experiences in similar soils in this region, with the main goal of defining the upper-bound and lower-bound estimates of the mass deposit extents and flow velocities of any debris flows originating from the steep unstable slopes. The bounds were kept deliberately large, given that no site-specific dataset were available. This reflected a high degree of uncertainty surrounding the choice of rheology parameters used in this assessment. The adopted rheological parameters are listed (Table 1). Figure 3 shows the relationship between the shear stress and shear rate characterized by the Herschel–Bulkley rheological model using three parameter sets.

Table  1.  Rheological parameters of debris flows

Case	Density/ (kg·m−3)	Herschel–Bulkley
		Yield stress/Pa	Consistency index/(Pa·sn)	Power index
LE (forward analyses)	1 524	160	6	0.67
BE (back and forward analyses)	1 695	900	21	0.56
HE (forward analyses)	1 695	1 640	21	0.56

 | Show Table

DownLoad: CSV

Figure  3.  Rheological relations of debris flows.

DownLoad: Full-Size Img PowerPoint

2.3.2   Pre-failure material source for back analysis

The back analysis simulation was conducted with the rationale that the measured present-day bathymetry was a post-failure surface affected by the relic mass movement events. Little evidence, however, linked the historic deposit to a source location upslope on the escarpment. Therefore, for this back analysis modeling effort, a hypothetical source location was selected toward the top of the escarpment within a deeply incised canyon (Fig. 4). The modeled source was a teardrop in plan view. The pre-failure bathymetry was constructed over the selected source region using the interpolation of the contour lines across the selected area. The pre-failure bathymetry within the boundary of a polygon delineated the edge of the failed soil masses. The bathymetry in the downslope zone was not changed. The differences between the reconstructed pre-failure surface and the present-day bathymetry resulted in a volume and thickness distribution of soil masses that was then simulated as participating in a relic debris flow event. The maximum thickness of the modeled source region was approximately 60 m with a total volume of 5 600 000 m³.

Figure  4.  Back analysis simulation: pre-failure bathymetry construction (a), and thickness distribution (b).

DownLoad: Full-Size Img PowerPoint

2.3.3   Pre-failure material source for forward analyses

This phase of work included the numerical modeling of debris flows originating from the slope stability profiles considered to be the most critical to the planned Manifold A, as identified through previous slope stability analyses. Figure 5 displays the selected locations of several slope profiles that were assessed during the slope stability analyses across the planned development area. The slope stability analyses demonstrated that Profile 2 and Profile 3, which were right next to each other, were the most critical slopes lying upslope of the Manifold A. They had identical factors of safety for the static condition and only a slight difference under pseudo-static conditions. The 2D and 3D dynamic finite element analyses further showed that the slopes would fail during or immediately following a seismic event, as enough elements experienced shear strains exceeding the point at which the strain-softening process was initiated. Subsequently, the unconstrained slope failures were likely to develop into a debris flow. Given the identical factors of safety for these two slope profiles, the implication was that both of the slopes and any other nearby regions with similar slope angles could simultaneously fail under these types of seismic conditions. Thus, the source facet for these forward prediction analyses was modeled such that it captured both slope profiles and other nearby areas with similar slope gradients. It was believed that this was a more realistic outcome of any seismic events compared with a small, localized failure.

Figure  5.  Slope stability profile: selected slopes for stability analyses (a), and example of slip surface (b).

DownLoad: Full-Size Img PowerPoint

The modeled sources for these forward analyses are presented in Fig. 6. The geometries of the source volume were derived qualitatively based on engineering judgments but informed by the previous results of 2D finite element dynamic slope stability analyses, along with the static and pseudo-static slope stability results. The slip surfaces predicted in the slope stability analyses were similar between the selected profiles, Profile 2 and Profile 3. The depth of the source volume was modeled to match the predicted slip surface at Profile 3. This slip surface was computed by examining the shear strains predicted in the 2D finite element dynamic slope stability analysis and assuming that all elements where these shear strains exceeded 5%. This was the strain level at which strain-softening was initiated. Naturally, the base of this region was the slip surface. This is illustrated as an example in Fig. 5b, where the red zone corresponds to elements where the predicted plastic shear strains exceeded 5%. The regional slope gradient within this canyon feature was used to map the potential slip surface and the extent of the source geometry. The maximum depth of the slip surface was usually accompanied by the highest slope gradient. Where the slope gradient was lower, the source thickness was reduced, in effect feathering the source toward the edge of the steep region. This formed an equivalent post-failure bathymetric surface and a credible 3D volume of the potential pre-failure mass source. The total volume of the material source modeled in this forward analysis was approximately 0.012 km³ (12 000 000 m³).

Figure  6.  Forward prediction simulation: post-failure bathymetry construction (a), and thickness distribution for the Manifold A (b) and Manifold B (c).

DownLoad: Full-Size Img PowerPoint

Based on the outcomes of the modeling and slope stability study as a whole, a similar analysis was performed with consideration of the slope stability profile, which was qualitatively considered to be the most critical to the proposed Manifold B and surrounding infrastructure. The pseudo-static slope stability analyses were not included for the slopes close to the Manifold B. This necessitated a slightly different approach in assessing the likely runout behavior of the debris flows originating upslope of the Manifold B and the potential risks they posed to downslope pipelines and fixed infrastructures. A deeply incised canyon feature was found to lie directly upslope of the Manifold B, and this was a likely source location for any debris flows developing from upslope failures in this region. The source volume for the debris flow runout simulations was derived based on engineering judgments and heavily informed by the slope stability results at Profile 3, along with the local slope gradients within this second canyon feature. The modeled material source is shown (Fig. 6). The total volume of the mass source modeled in the forward analysis was approximately 1 150 000 m³, which was 10 % of the size of the examined debris flow originating upslope of the Manifold A. Note that, these debris flow runout events were modeled by adopting the conservative assumption that, following an initial slope failure, the entire soil mass over the slip surface disintegrated and developed into debris flow as it moved downslope.

3. Results and analyses

3.1 Back analysis runout modeling results

Figure 7 illustrates the resulting mass deposit obtained when the simulation was run using the BE rheological parameters listed in Table 1 and the source presented in Fig. 4. The figure illustrates the flow path and thickness of the deposit once it came to rest, as well as the maximum velocities encountered over the course of the simulation. After the initial phase of the runout in which the flow was characterized by the material moving off the canyon flanks into the thalweg, the debris flows tended toward a south-east direction down the slope. Flow deposits of up to 2.5 m thick were found within the runout lobe and along the flow track edges at the base of the steep slope. The flow velocity peaked at the base of the steep slope, at approximately 50 m/s. The computed runout extended slightly beyond the extent of the mapped MTD. This might have been because the MTD likely continued significantly beyond its mapped extent with a reduced thickness that could not be fully resolved in the 3D seismic data. Whilst increasing the yield stress of the soil would increase its resistance to the flow and better match the extent of the observed debris mass deposit, it was difficult to justify this without site-specific rheology testing or better resolving the extents of the mapped MTD through the acquisition of improved shallow geophysics.

Figure  7.  Back analysis flow path, thickness, and velocity following the runout event.

DownLoad: Full-Size Img PowerPoint

3.2 Forward analyses runout modeling results

3.2.1   Prediction results for the Manifold A

Figure 8 illustrates the resulting mass deposit obtained when running the analyses using the LE, BE, and HE rheological parameters listed in Table 1, with the source presented in Fig. 6. The figure illustrates the flow path and the thickness following the event, as well as the maximum velocities encountered over the course of the simulation. There was a moderate sensitivity to the LE, BE, and HE rheology parameters in each simulation. The extent of the flow deposit varied somewhat depending on the rheology parameter used, with the LE and BE cases extending farther down the canyon compared with the HE case. The runout distances, which were measured from the source to the nose of the debris flow, ranged from 5 km with HE rheology to 15 km when LE rheology parameters were adopted.

DownLoad: Full-Size Img PowerPoint

Figure  8.  Forward predictions of flow path, thickness, and velocity following the runout event for the Manifold A: Low Estimate (LE) case (a), Best Estimate (BE) case (b), High Estimate (HE) case (c), and flow extent summary (d).

DownLoad: Full-Size Img PowerPoint

Table 2 summarizes the approximate area covered by the simulated debris flow and minimum distance to the Manifolds A and C for each of the LE, BE, and HE rheology variations. It was clear from the modeled results that the debris flow originating from the source location might pose a potential risk to the infrastructure near the proposed Manifold A, as the predicted flow paths in both the LE and BE cases crossed the proposed pipeline alignment about 1 km to the west of the manifold location. The HE case stopped short of reaching the proposed structure with a runout distance of about 5 km. While the simulations had different mass deposit extents, the general directions and flow velocities of each flow were very similar. All three cases considered in this study tended in the south−east direction down the slope following an initial phase of flow characterized by the movement of source materials off the canyon flanks into the canyon thalweg.

Table  2.  Simulated area of mass deposit and minimal distances to the Manifolds A, B, and C

Rheology variation	Simulated area/km2	Minimum distance/km	Manifold type
Low Estimate (LE)	14.3	0.7
Best Estimate (BE)	4.76	1.4	Manifold A
High Estimate (HE)	2.9	2.7
Low Estimate (LE)	3.97	4.5
Best Estimate (BE)	1.47	6.9	Manifold B
High Estimate (HE)	0.89	7.8
Low Estimate (LE)	14.3	2.9
Best Estimate (BE)	4.76	4.7	Manifold C
High Estimate (HE)	2.9	6.2

 | Show Table

DownLoad: CSV

The deposits within the runout lobe and along the flow track edges also varied in terms of the thickness and maximum velocities encountered during the course of simulation depending on the rheological parameters used. The maximum flow thickness developed in the HE case was about 2.3 m, with a maximum velocity of approximately 50 m/s located at the base of the source region as the flow moved through the canyon thalweg. In the LE case, relatively similar velocities were observed as the flow entered the canyon base and traversed the planned infrastructure. In the BE case, a velocity of approximately 15 m/s was observed as the flow traversed the planned infrastructure. The debris flow finally came to rest with the nose extending approximately 370 m past the pipeline alignment. Deposits up to approximately 20 m thick in the HE flow and 18 m thick in the BE flow were found to accrue in local topographic lows on the canyon floor and come to rest in a deposit fan. The thickest deposits were confined to the source area, which appeared to be caused by the ponding effect within this region, in which the source materials failed to gain sufficient energy during the initial acceleration to fully evacuate from the source region. No such deposits were observed in the LE case. As expected, the scenario adopting HE rheology parameters displayed shorter runout distances but thicker deposits closer to the source of the flow when compared to the scenarios adopting LE and BE parameters. However, the differences between the HE and BE cases were much smaller in comparison to the LE and BE cases.

3.2.2   Prediction results for the Manifold B

Figure 9 illustrates the mass deposit results obtained when running the analyses using the LE, BE, and HE rheological parameters listed in Table 1 and the material source volume illustrated in Fig. 6. The figure illustrates the flow path and the deposit thickness following the flow event, as well as the maximum velocities encountered over the course of the simulation. As seen with the Manifold B simulations, each simulation had moderate sensitivity to rheological parameter variations, with the flow deposit for the LE case extending farther down the canyon compared with the BE and HE cases. The runout distances from the source location to the head of the flow were 6.5 km, 3.0 km, and 2.0 km in each of the LE, BE, and HE cases. Table 2 summarizes the approximate area covered by the simulated flows and the approximate minimum distance to the Manifold B for the LE, BE, and HE rheology variations. Although the simulations had slightly different deposit extents, the general directions and flow velocities of each flow were again quite similar. For all three cases considered in this study, after the initial runout, the flows tended toward the south−east direction down the slope.

Figure  9.  Forward predictions of flow path, thickness, and velocity following the runout event for the Manifold B: Low Estimate (LE) case (a), Best Estimate (BE) case (b), High Estimate (HE) case (c), and flow extent summary (d).

DownLoad: Full-Size Img PowerPoint

The flow deposits within the runout lobe and along the flow track edges also varied in thickness and maximum velocities depending on the rheological parameters used. For the debris flow simulated using the LE rheology parameters, velocities of approximately 0.5 m/s were observed toward the front tip of the flow. For the BE and HE cases, velocities observed toward the front tip of the flow were higher than in the LE case, with the BE case velocities of 1.3 m/s and 2 m/s for the HE case. When the HE rheology parameters were used, flow deposit thicknesses of up to 3 m were found to occur in the local topographic lows on the canyon floor and as the flow came to rest in a deposit fan. As expected, the scenario adopting the HE rheology parameters displayed a shorter runout distance and thicker deposits than the scenarios adopting the LE and BE rheological parameters. In all the scenarios, the flows were predicted not to pass over any of the planned infrastructures surrounding the Manifold B. Thus, the simulations suggested that the planned infrastructures would face little risks from future potential debris flows with similar magnitudes originating from the upslope incised canyon.

4. Discussion

4.1 Limitations of the modeling results

Modeling the runout of submarine debris flow is a challenge because the hydrodynamic processes of such a real-world flow are quite intricate (Mohrig et al., 1998; Marr et al., 2001; Toniolo et al., 2004; Ilstad et al., 2004a, b; Elverhøi et al., 2005, 2010; Zakeri et al., 2008; Haza et al., 2013; Du et al., 2022). In this study, the modeling results of submarine debris flow were based on a numerical code developed in-house by Fugro (Spinewine et al., 2011). To simplify the problem, a series of assumptions was made, and some complex physical mechanisms were not explicitly considered in this model. As a result, simulation results from this model may suffer from some significant limitations.

Faced with the problem of assessing potential runout of a submarine debris flow, it is commonly assumed that the initial strength of a debris flow is equal to the fully remolded material strength. However, the values of the material strength required to match the runout of historical failures are often an order of magnitude lower than the fully remolded debris strength (Imran et al., 2001; De Blasio et al., 2004; Niedoroda et al., 2006; Das, 2012; Qian et al., 2020). Hence, the best-fit model parameters are often poorly related to the soil parameters obtained through conventional geotechnical tests (Ingarfield et al., 2016). In addition, the shear strength of debris flow is assumed to remain constant during the runout process. The effects of phenomena, such as gradual changes of material strength caused by either water or soil entrainment (De Blasio et al., 2005; Gauer et al., 2005; Kim et al., 2019), are not considered explicitly, and they are incorporated implicitly through model calibration exercises of historical failures. Based on this setting of material strength in the numerical model, the runout distances of past events are usually well calibrated (Imran et al., 2001; Das, 2012; Qian et al., 2020). However, the predicted velocities of debris flows will still experience some kinds of uncertainties and might deviate significantly from the reality.

The Herschel–Bulkley fluid flowing down a slope is characterized by a basal shear layer overlain by a rigid plug without an internal shear rate, and the model assumes a no-slip and rigid bed condition. Material exchanges with the seafloor in terms of further erosion and deposition outside of the initial source area were neglected in this study. Because deposition was not considered explicitly, termination of the debris flow occurred only when the basal shear stress fell below the prescribed yield stress. The ignorance of material exchanges with the seabed typically in the form of substrate erosion and debris mass deposition will also lead to significant uncertainties in reliably modelling the runout of a debris flow (Mohrig et al., 1999; Toniolo et al., 2004; Nugraha et al., 2022; Du et al., 2022). For example, significant erosions of substrate will bulk up the total volume of propagating debris material, whereas deposition processes will reduce the amount of moving material. As a result, the total volume of a debris flow and associated rheological parameters will be significantly modified, which leads to further modifications of the runout distances and accompanying velocities.

In reality, a real-world submarine debris flow is thought to consist of two distinct layers, which is a high-density laminar layer overlain by a much lower-density turbulent layer of turbidity current (Mohrig et al., 1998; Marr et al., 2001; Toniolo et al., 2004; Ilstad et al., 2004a, b; Elverhøi et al., 2005, 2010; Zakeri et al., 2008; Haza et al., 2013; Du et al., 2022). The Reynolds number is used to help inform if the debris flow as a whole is in laminar, turbulent, or transitional flow regimes, although the boundaries between these three regimes are not well defined (Elverhøi et al., 2005, 2010; Zakeri et al., 2008). For the large, deep-rooted debris flows reaching high velocities, the Reynolds number typically will show that the flow is in laminar or transitional flow regime, which is consistent with the two-layer understanding. In this scenario, shallow and soft layers of soil at the surface will tend to mix up with the ambient seawater and therefore become fully turbulent. This type of turbulent flow would act quite differently from the modeled flows in this study. As a consequence, this turbulent weak layer at the top would tend to trail behind the leading edge of a debris flow and not actually participate in its overall dynamics (Marr et al., 2001; Toniolo et al., 2004; Ilstad et al., 2004a, b). In contrast, the basal layer of the flow head, consisting of denser and stiffer material from the failed masses, would be likely to remain quasi-laminar, in agreement with the modeling assumptions. The modeling in this study considered only the underlying high-density laminar layer, ignoring the material exchanges with the overriding turbidity current. In this respect, the modeling could be thought of as somewhat conservative. Both mixing up with and inertial resistance from ambient seawater tended to decelerate the debris-flow velocity (Gauer et al., 2005; Kim et al., 2019; Qian et al., 2020). This deceleration effect was the most pronounced in the region of highest debris-flow velocity, especially in frontal heads where stagnation pressure from ambient seawater was the largest. Therefore, this modeling framework was prone to overestimating the frontal debris flow velocity. Hence, it must be kept in mind that the values and maps of the maximal flow velocity reported in this study were likely to be overestimated. In contrast, the velocities of main body and the trailing edges of debris flow were much less affected by the phenomena and were expected to be better simulated in this model.

In addition, this numerical model did not consider the possibility of hydroplaning of submarine debris flow, which could lead to much higher velocities being experienced, because a thin layer of water underlying the head of a debris flow acted to lubricate the movement of the soil mass (Mohrig et al., 1998; Marr et al., 2001; Toniolo et al., 2004; Ilstad et al., 2004a, b; Elverhøi et al., 2005, 2010; Zakeri et al., 2008; Haza et al., 2013; Du et al., 2022). The mechanism was not sufficiently understood at the time of the study to incorporate it into the numerical code. However, ignoring hydroplaning balanced the conservative assumption of ignoring the resistances from water pressures imposed on the leading edge of the debris flow to some extent.

Finally, the limitations of the in-house developed model (Spinewine et al., 2011) pointed out as above still remain to be addressed in future work, and there is a significant room to advance the state-of-art runout modeling technology of large-scale submarine debris-flow scenarios from the perspective of real-world offshore engineering practices. Despite the limitations, the modeling results presented in this study can still serve as a preliminary reference to evaluate submarine debris-flow hazard risks to planned subsea pipeline systems in this development region, and provide scientific support to the design of a more resilient pipeline alignment in consideration of future debris-flow hazard risks.

4.2 Optimizations of the pipeline route

The runout analyses of submarine debris flow discussed in this study showed the possibility of interactions between planned subsea pipelines and potential future debris flows in the development site. Several avenues could be explored to demonstrate a reduced hazard risk to these pipelines from future mass transport events. One of these avenues is to alter the planned alignment of these pipelines to avoid the likely flow path of a debris flow. Avoiding the worst of these debris flows, that is, the debris flows simulated using LE rheology parameters, is likely impractical given the extensive runout predicted in this case. It is likely possible to reroute these pipelines to avoid the potential runout path of the simulated debris flow with BE rheology parameters and to minimize exposure to any other flows of this magnitude originating from the upslope escarpment. Herein, two possible amended pipeline alignments are recommended as follows.

If it is feasible that the Manifold A be relocated farther downslope, then it could be brought in line with the Manifolds B and C (Fig. 10a). This alternative would provide the most benefit in mitigating the hazard risks posed by any potential debris flows. Pipelines crossing the likely path of any runout originating from the deeply incised canyon would be located significantly farther downslope than the toe of mass deposits simulated with the BE rheology parameters and would continue to the Manifold B at a location farther away from the base of escarpment when compared to the originally planned alignment.

Figure  10.  Pipeline route optimization: alternative 1 (a), and alternative 2 (b). HE: High Estimate; BE: Best Estimate; LE: Low Estimate; ILT: in-line-T.

DownLoad: Full-Size Img PowerPoint

Alternatively, if the Manifold A cannot be relocated, for example, if wells are located nearby, then an alternative strategy could be pursued. The option would leave the Manifold A at its original location, but would reroute the export pipeline farther downslope to intersect the Manifold C (Fig. 10b). From this point of the Manifold C, the pipeline initially would travel directly toward the Manifold B before branching off and tracking back upslope to intersect the Manifold A. Although the Manifold A is left closer to the submarine escarpment than with the previous alternative, it is to the east of a small topographic high that would help deflect any debris flow away from the manifold should debris flow sourced from the deeply incised canyon occur. This topographic feature, however, would offer no protection from debris flows originating farther east. In addition, an in-line-T (ILT) structure would be required such that infield flowlines can be used between the ILT and Manifolds A and B, while the flowline would be maintained from the ILT to the Manifold C.

5. Conclusions

In this case study, a quantitative evaluation of the submarine debris-flow hazard risks to planned pipeline systems in a new development site in the Qiongdongnan Basin, South China Sea, was performed using a numerical program developed in-house by Fugro. Runout simulations of the debris flow were performed adopting the LE, BE, and HE rheology parameters to bracket different extremes of potential flow behaviors at the development site. The modeling results showed that a significant runout distance of a potential debris flow could be expected. In particular, when using either the LE or BE rheology parameters, the debris flows could interact with the planned pipeline systems around the Manifold A, and these pipelines would likely not be survivable under this impact. It is worth noting that these analyses adopted several simplifying assumptions out of necessity and deliberately erred on the conservative side. Even so, the results presented in this case study can still serve as a reference to support the design of a more resilient pipeline alignment in consideration of future debris-flow hazard risks. Based on the debris flow modeling analyses, it is recommended that consideration be given to two possible alternative pipeline alignments that will reroute the proposed pipelines much farther downslope and out of the runout path of the best-estimate flow.