Numerical Study on Interaction between Submarine Landslides and a Monopile Using CFD Techniques

: Offshore installations with pile foundations in shallow water are vulnerable to submarine landslides, which cause serious damage to engineering facilities, loss of life, and loss of money. Due to a shortage of real observation data and the difﬁculty of reproduction, we lack insight into the interaction behavior between submarine landslides and monopiles. This study capitalized on ANSYS Fluent 20.0 to develop a three-dimensional biphasic (water and slurry) numerical model. This CFD model was used to analyze the interaction between a monopile and submarine landslides at different ﬂow heights. The velocities of submarine landslides were from low to high values. Two modes of interactional forces acting on the monopile are proposed, which are (i) interaction force with peak value and (ii) interaction force without peak value. The inﬂuence of ﬂow height and velocity on interaction forces was investigated. Results show that the effect of the ﬂow heights on the interaction force is signiﬁcant at low velocity stage, while the peak force representing a hazard level of the pile was non-negligible under high ﬂow velocity and low ﬂow height conditions, which should be considered in a future study. The related mechanisms are revealed with a hybrid model considering different components of the force.


Introduction
With the increasing exploitation of offshore resources, many offshore installations have been put into use, such as fixed platforms, submarine pipelines, and other subsea facilities [1]. These installations may be exposed to the devastating threats of submarine landslides, which are a common and frequent marine geohazard [2,3]. Submarine landslides may be triggered by many factors [4], such as volcanoes and earthquakes, waves, dissociation of hydrate methane, and so on. Among them, storm surge is a very possible trigger mechanism in shallow water area [5][6][7]. Under the influence of wave dynamic loading of the storm surge, seabed instability may happen in this area [8]. When the seabed becomes unstable, submarine landslides can gradually disintegrate into fluidized debris flows [9], with velocity up to 20 m/s [10][11][12][13], which have the potential to severely damage structures along the paths. This type of hazard usually causes injury or loss of life and imposes significant additional project costs for mitigating their effects [14][15][16]. In shallow water areas, piles are widely used as foundation for engineering structures [1,17,18]. To ensure the safety of design and operation of offshore structures with pile foundations, estimating the potential interaction force between submarine debris flows and a monopile is of great research significance and practical value for marine engineers.
Recently, a large number of studies on submarine landslide-structure interaction have been carried out. Among them, the interaction between the submarine landslides and pipelines have particularly attracted researchers' attention. In these studies, the major focus is the debris flow stage of submarine landslides, which has the potential for catastrophic disasters to seafloor structures. For the debris flow stage, fluid mechanics principles seems to be more appropriate [10,19,20]. Based on these theoretical principles, flume tests were performed to estimate the interaction force acting on pipelines, which lays the foundation in this research field [21]. Then CFD analyses [22] were conducted to complement the flume tests results. The data of CFD calculations combined with the results from the flume tests are meant to propose the relationship between the drag force coefficient and Reynolds number. Next, some researchers extended the previous works to study the normal and longitudinal interaction forces at different impact angles of submarine debris flows [23][24][25]. Afterwards, more factors on the drag force, such as suspension height [26][27][28], temperature of submarine debris flows [29], shape of the pipeline [30][31][32], strain-softening effects on soil strength [33], were also considered. In addition, some studies [34,35] used the material point method to investigate impact forces on pipelines, while an explanation of the mechanism of peak force was provided by [36,37].
It should be noted that previous studies mainly focused on the interaction force acting on pipelines, which is different from that acting on a pile. Feng et al. [38] launched a small number of flume tests to investigate the interaction force acting on a monopile at moderate flow velocities (0.38-0.94 m/s) through preliminary research. However, due to the limitation of the Feng's experiment, tests in high-speed and low-speed movement of submarine landslides were not considered. Moreover, in offshore environments, flow heights of submarine debris flow have a remarkable influence on the force acting on a monopile, and the effect of flow heights has not been discussed in the flume tests. The understanding of the interaction processes and underlying mechanism is still inadequate due to insufficient data.
In this paper, a three-dimensional biphasic (slurry and water) numerical model using ANSYS Fluent 20.0 was developed. This model, based on the Herschel-Bulkley rheology, was used to investigate the interaction between submarine debris flows and a monopile. Firstly, the numerical model was verified by the results of a flume test. Then, a large number of simulations were carried out and two models of interaction force are summarized. The effects of flow velocities and flow heights on interaction forces on the monopile are also discussed. A qualitative variation law of the interaction force is proposed and the mechanism underlying is revealed through a hybrid model considering different components of the force. In addition, the difference between the peak and the stable drag coefficients is discussed to estimate the influence of the peak value, which is highly harmful for a monopile.

Interaction Model
In shallow water areas, after the seabed failure is triggered by large waves, submarine landslides may flow downslope. In this process, the submarine landslide body may evolve into submarine debris flows, which have greater influence and larger potential damage for marine structures. The present research mainly focus on the submarine debris flow stage, where the landslide body material may be considered as a fully remolded viscous fluid [9,39]. The interaction force between submarine debris flows and the monopile plays a vital role in the safety of offshore installations. In this study, due to the complexity of the interaction force in the underwater environment, the simulation just focuses on the force acting on the monopile. The monopile was fixed on the seabed slope whose surface was a rough wall boundary. Moreover, it is assumed that the monopile does not deform when interacting with submarine debris flows, as shown in Figure 1. This is therefore an upper bound of the loading scenario; that is, the estimation of the interaction force is relatively safe.
The interaction force of submarine debris flows acting on a monopile is divided into three components: drag force (normal), friction force (axis), and transverse force (lateral). The primary concern in this study is the drag force, which represents the most detrimental load. For the fluidized debris flow stage, the drag force exerted by the submarine debris flows and water on the monopile was calculated in the simulations using the fluid dynamics method, which is simpler and more popular [40],. Therefore, the drag force coefficient, C D , is defined in the same way as that in Zakeri et al. [24].
where C D , ρ and D are the drag force coefficient, the density of a submarine debris flows, and the diameter of a monopile, respectively. A is the normal projected frontal areas of a monopile. H is the flow height of submarine debris flows. U ∞ is the approaching velocity of the submarine debris flows head, which was taken as the flow velocity of submarine debris flows. This is because when hydroplaning occurs, the submarine debris flows head advances as a plug flow with a uniform velocity profile [24]. F D is the measured total drag force. For the Reynolds number of the submarine debris flows (Re non-Newtonian ), the definition, with τ as the shear stress of the submarine debris flows, is given by Equation (2).

General
Currently, many numerical methods, which can deal with large deformation problems, are available for simulating the post-failure runout of submarine landslides. Such numerical methods include: the coupled Eulerian-Lagrangian method (CEL) [41], mesh-based large-deformation finite element analysis (LDFE) [42], meshless material point method (MPM) [35], smoothed particle hydrodynamics (SPH) [43,44], computational fluid dynamics with the discrete element method (CFD-DEM) [45], and so on. Among them, computational fluid dynamics (CFD) is one of the most widely used numerical methods [46][47][48][49][50]. The validity of a numerical method for calculating the impact of submarine debris flows on pipelines has been successfully demonstrated by using CFD [22][23][24]36,51]. The CFD software ANSYS Fluent 20.0 (Ansys student 2020 R2 (Built-in license valid until 7/31/21).) was employed in this study, which is a general-purpose CFD software that includes solver based on the finite volume method (FVM) for unstructured meshes.
The submarine debris flows constitutes an incompressible biphasic (slurry and water) flow regime. As such, the Eulerian multiphase model is used to simulate the constant-speed fluidized submarine debris flows interacting with a monopile. In the numerical model, conservation equations of mass and Navier-Stokes equations method are discretized using the FVM. They are solved separately for each phase. The interphase coupling is achieved through the shared pressure and interphase exchange coefficients. A general descriptions of the theory and associated formulations regarding Eulerian multiphase model used to analyze the result can be found in the ANSYS Fluent Theory Guide [52].

CFD Model Setup
Feng's experimental data [38] were complemented using a similar CFD model by simulating mudflow entering the domain full of water at various flow velocities and heights. Parameters of submarine debris flows in the numerical model was strictly within the applicable range, which does not affect the fundamental understanding of the qualitative variation law and mechanism underlying [36,53]. The numerical model used to simulate the submarine debris flows-pile interaction is shown in Figure 2. The pile diameter (D) used in the simulations was 32 mm. The domain is 11D wide, 31D long, and the maximum height is 31D at the shallowest point (i.e., at the outlet). The slope angle of seabed is set as 4 degree. Figure 2b shows the pile is 12.5D high, whose top surface is 8D away from the top boundary of the domain. The monopile was more than 5D away from all the boundaries of the domain, which was found through sensitivity studies to be sufficient to eliminate boundary effects.
Both the monopile and seabed surfaces were modeled as rough no slip boundary conditions with the equivalent sand-grain roughness, k s , of 1.6 × 10 −6 m [54] and 12.5 × 10 −6 m [38], respectively. The submarine debris flows of different flow heights (H) entered from the left side of the domain, with a constant velocity parallel to the seabed, as depicted in Figure 2c. The top and both sides of the domain as well as the wall above the inlet were modeled as free-slip wall boundary condition, which reduces the computation time required to compute the water flow field on these surfaces without affecting the submarine debris flows [24]. The outlet was modeled as outflow boundary condition. The domain is discretized into unstructured meshes using the ANSYS Meshing module, as shown in Figure 3. The average maximum mesh size was 3.5 × 10 −2 m. Since submarine debris flow-pile interaction is of main interest for the study of fluid-structure interaction, meshes around the monopile were refined to a maximum mesh size of 0.01 m. The inflation boundary of five layers with growth rate of 1.2 was set at the pile edge. The maximum mesh size on the monopile surface was 2.5 × 10 −3 m and gradually increased to 0.01 m within the 3D outward radial distance. Through the mesh sensitivity analysis, a mesh of 390,000 magnitude was confirmed for the calculation to balance the computational efficiency and accuracy, as shown in Figure 4.   Table 1 and standard properties of water were used in the simulations, respectively. Turbulence generated in the two-phase flow was simulated using the κ-ω model. The difference between two-phase flow velocities produces phase interaction at the contact surface. The phase interaction was set as drag coefficient. The total time of simulation was 3-10 s (related to mud velocity), and the time step was 0.0005 s. The transient formulation and spatial discretization were taken as the second order form. Note: the Herschel-Bulkley model is expressed with τ = τ y + Kγ n . τ is the shear stress, τ y is the yield stress, K is the consistency index, n is the flow behavior index and γ is the shear rate.

Comparison and Validation
To verify the validity and reliability of the numerical calculations, various operating conditions of Feng's study [38] were simulated. Figure 5 shows the comparison between the current calculations and previous results. Generally speaking, the current drag force coefficients are quite similar to the results obtained in the flume tests. The simulation yielded higher drag force coefficients than that in the experiment at low Re non-Newtonian . The reason for this discrepancy is that the mudflow slightly slip on some parts of the pile in the flume experiment, whereas no slippage is allowed in the numerical simulation. As for high Re non-Newtonian , the numerical results are lower than the experimental data. This may be due to the large discreteness of experimental results [30].

Load Characteristic of the Pile
The simulations were extended to include the full range of flow velocities. The effect of different flow heights on the drag force was systematically studied. The rheological parameters and the seabed slope utilized in all of the simulations were the same as the ones used in the flume experiments. The ranges of submarine debris flow parameters used in the simulations are shown in Table 2. A total of 85 CFD simulations were conducted.  Figure 6 shows a typical loading response imposed by submarine debris flows on the monopile. For low U ∞ and H, no peak drag force is observed. When U ∞ and H increase, the peak value appears. Two typical models can be summarized from all the simulations conducted, which are model 1 (high U ∞ , H) with a peak drag force value and model 2 (low U ∞ , H) without a peak drag force value. We take submarine debris flows with U ∞ = 0.3 m/s, H = D (see curve model 2 in Figure 6) and U ∞ = 1.25 m/s, H = 5D (see curve model 1 in Figure 6) as two examples to illustrate both models. For model 1, the drag force increases with time sharply up to point A 1 and then decreases to a negative value. Then the value quickly rises to a peak at point C 1 , and afterwards gradually declines until it reaches a stable value at point D 1 . For model 2, the stage before point B 2 is similar to model 1. After point B 2 , the drag force gradually increases to a stable value at point C 2 . The formation mechanism of the models is explained as follows. At point A 1 and A 2 , the landslide body with a certain initial velocity enters the domain from the inlet and pushes the seawater in front so that the originally stationary seawater begins to accelerate to a higher velocity and first impacts the monopile, causing a sharp increase of the force [31]. As the landslide moves forward, the left end and right ends of the pile can be divided into upstream and downstream zones. For model 1, before the debris flows touch the monopile, strong vortices form upstream of the pile, resulting in a negative pressure zone. This slows down the growth of the drag force before point B 1 , as shown in Figure 7a. After the landslide body gradually comes into contact with the monopile, the pressure difference between the positive pressure zone and the negative pressure zone increases and the drag force gradually increases to its peak value, as shown in Figure 7b. As more of the landslide body engulfs the monopile, the negative pressure zone becomes narrower. The drag force decreases to the stable value, as shown in Figure 7c. For model 2, since the velocity and flow height are lower, it takes longer for the landslide body to meet with the monopile at point B 2 and the drag force is smaller than that of model 1, as shown in Figure 7d. It can be seen from Figure 7e that drag force increases to the stable immediately after point B2 because the pressure difference between the upstream and downstream zones is negligible. Generally, in addition to the velocity of submarine debris flows, the flow height influences the drag force acting on the monopile. Therefore, it is necessary to perform a detailed analysis to reveal the underlying law and mechanism.

Influence of the Submarine Debris Flow Heights
Considering the different flow heights (H), the computed drag force was analyzed utilizing the same method as described in Section 2. Figure 8 shows the relationship between C D and Re non-Newtonian at different H. The peak drag force coefficient (C D-P ) and the stable drag force coefficient (C D-S ) decrease with Re. A critical Reynolds number (Re c = 165.5, U ∞ = 1 m/s) is observed. The overall trend at low Reynolds number stage (Re non-Newtonian < Re c ) and high Reynolds number stage (Re non-Newtonian > Re c ) are different.  Figure 9a shows that at the lower Reynolds number stage (Re non-Newtonian < Re c ), C D-P increases with H and reaches the maximum value at H = 9D. C D-P remains almost constant and the gradient of C D-P equals to zero, when Re non-Newtonian = Re c . When Re non-Newtonian > Re c , C D-P decreases with H to a minimum at H = 9D. Referring to Figure 9b, C D-S gives a similar trend as that of C D-P at the lower Reynolds number stage (Re non-Newtonian < Re c ). However, C D-S for the full range of heights remains the same after Re non-Newtonian is greater than Re c . In addition, the gradient of C D along with H when Re non-Newtonian < Re c is greater than that of higher Re non-Newtonian , and when Re non-Newtonian > Re c , the gradient of C D is much smaller than that of the lower Reynolds number range. In short, the effect of H on C D is significant at the lower Reynolds number stage (Re non-Newtonian < Re c ). According to the formation mechanism [35], the drag force is composed of the differential pressure inertia force, the viscous strength force and the static pressure force. The viscous strength force is closely related to the roughness of the monopile surface and the debris flows movement within the boundary layer. However, since the same flow model is used in this study and the roughness of the monopile is consistent in all simulations, the contribution of the viscous strength force is not considered. For low Re non-Newtonian conditions, the static pressure force closely related to the flow height (H) is dominant. Therefore, the drag force increases with H. For high Re non-Newtonian conditions, the static pressure force is negligible, hence the differential pressure inertia force, which is due to the separated region of the monopile, is considered [30]. Figure 10 compares the mud volume fraction and pressure distribution of the same submarine debris flow with the same high impacting velocity at the time of the peak drag force. At this moment, the shape of the downstream zone for H = D is significantly different from those for H = 5D and 9D, resulting in the different separated region. The separated region has a direct influence on the negative pressure zone, which leads to the variation of the differential pressure inertia force.  To understand the quantitative variation of the drag coefficients, the maximum rate of change of C D (β) with different H is obtained, as shown in Figure 11. When Re non-Newtonian < Re c , β of C D-P is greater than that of C D-S . The maximums of C D-S and C D-P are 54% and 91%, respectively. Basically, the influence of H on C D-P is greater than C D-S , and the difference between C D-P and C D-S should be elaborated.  Figure 12 shows the variation of the residual of drag force coefficients (R) with respect to Re non-Newtonian . The overall trend shows that R decreases with Re non-Newtonian until it reaches a critical value at Re C = 165.5, and then begins to increase. For H = D, when Re non-Newtonian < 20, R is close to 1, indicating that the influence of H is negligible, which is consistent with the law of model 2. The focus in this study falls within the low Re non-Newtonian range (Re non-Newtonian < Re c ) because of the significant impact on C D . When Re non-Newtonian < 100, a greater H of submarine debris flow corresponds to a greater R, which means a greater gap between the peak drag force and the stable drag force. However, considering the flow height (H), the variation law of R in the middle Re range (100 < Re non-Newtonian < 250) is opposite that in the lower Re range (Re non-Newtonian < 100). For the middle Re, the minimum (R = 0.638) appears at H = D, which the means peak drag force is more significant with higher Re non-Newtonian and lower H. In a previous study [38], due to the limitation of the tests, the stable drag force, which lasts a relatively long period, was selected to analyze the drag force. However, the present study shows the peak drag force is significant and highly harmful for the pile.

Conclusions
In this paper, the interaction force between submarine debris flows and a monopile was systematically analyzed with computational fluid dynamics (CFD). To complement flume experiment data, flow velocities of the submarine debris flows were varied in the range of 0.1-2.5 m/s. Moreover, the effect of the flow height was carefully investigated. The main conclusions are as follows.
(1) The numerical results agreed well with the experimental results. Therefore, such a CFD method may be utilized to investigate the interaction between the submarine debris flows and a monopile with similar settings. (2) The typical load characteristics of the drag force acting on the monopile were analyzed.
Two modes of the drag force were proposed. The flow velocity and the flow height of submarine debris flows were found to have an influence on the mode of the drag force. The underlying mechanism was revealed using a hybrid model combining the effects resulting from submarine debris flows differential pressure inertia force, viscous strength force and static pressure force. (3) For low Reynolds number (Re non-Newtonian < Re c ), the flow height of the submarine debris flows (H) has a non-negligible influence on C D , which should be a major consideration in future studies. Moreover, the influence of H on C D-P and C D-S is different. The maximum rates of change of C D-S and C D-P with respect to different H can reach up to 54% and 91%, respectively. (4) The gap between C D-P and C D-S changes with Re and H. The maximum gap appears at Re c = 165.5 and H = D in this study. That is, in engineering design, the peak value, which represents the hazard level of the monopile, may also need to be considered.
The discussion on the interaction force in this study is preliminary, and some issues need to be addressed further. Due to the complexity of the marine geological environment, a simplified interaction model is adopted to estimate the interaction force, which is relatively conservative and limited. Some factors, such as the deformation of the monopile, the interaction between the seabed sediment and the monopile, as well as the applicability of parameters, have not been considered. To obtain a more comprehensive interaction model, additional physical model tests and large-scale simulations need to be performed in the future.