Rheological considerations for the modelling of submarine sliding at Rockall Bank, NE Atlantic Ocean

Recent scientific research indicates that the Rockall Bank Slide Complex in the NE Atlantic Ocean has formed as the result of repetitive slope failures that can be distinguished in at least three major phases. These sliding episodes took place during and before the Last Glacial Maximum. This work attempts the modelling of each sliding episode with the incorporation of the landslide’s rheological properties. The objective is to study the landslide kinematics and final deposition of each episode under a rheological framework that comes in agreement with the field observations. To do so in the present work, we use different types of rheological models to compute the total retarding stress and simulate submarine failure. The Bingham rheology and the frictional rheology are used to model the flow behavior. The scope of this approach is to understand the effect of the two classical laws in landslide kinematics. A rheological model that combines the two regimes is also used. To account for the hydrodynamic drag,...


I. INTRODUCTION A. Rheological regimes for the description of submarine motion
The post-failure stage of submarine landslides can be adequately described using the fluid mechanics principles and is often considered as the most significant for tsunami generation. 1,2Above all, the sediment concentration and water content of the landslide play an important role on the rheological behavior.For example, landslides with water content greater than 50% behave as Newtonian fluids (classified as streamflows), whereas landslides with less water content exhibit a visco-plastic rheological behavior (classified as slurry and granular flows). 3epending on the selection of the input parameters for the calculation of the shear retarding stress, the flow properties may vary corresponding to different fluid types.Three main types of rheological models have been used in the past to describe the dynamics of submarine landslides (considered as non-Newtonian fluids) at the post-failure stage: 1,4 (1) viscous models, (2) frictional models, and (3) viscoplastic models.De Blasio et al. 5 referred to the latter two as the main flow models for the simulation of submarine motion.
The appropriate rheological regimes that can be used in the simulations may vary between different environmental settings depending on various factors.Among those are the sediment composition, the particle size distribution and shape, and the pore pressure. 6,7For example, three types of fluids have been commonly used to represent the rheology of mudflows or muddy debris flows (silt-rich and clay-rich flows): a Bingham fluid, a Herschel-Bulkley fluid, and a bilinear fluid. 1 On the other hand, terrestrial events and sand-rich flows exhibit a granular behavior where there is a lot of grain-to-grain interaction.These types of flows are better described by frictional flow models. 5,7,8ollowing the description of the rheological regimes as found in Ref. 2, the Bingham rheology yields where T is the magnitude of the retarding stress of the flow, T 0 is the yield strength, representing a critical threshold beyond which the body begins to flow, µ is the dynamic viscosity, and γ is the shear rate (or also strain rate).For zero viscosity, the second term is neglected and the flow motion becomes subject to a constant retarding stress. 9he Herschel-Bulkley rheology yields where K is a consistency factor and n is the flow index.For n < 1, the state of mixture describes a shear-thinning fluid (or else a pseudoplastic fluid), a fact which indicates that the viscosity decreases as the shear rate increases.On the other hand, for a shear-thickening fluid (or else dilatant fluid) where n > 1, the viscosity behaves in the opposite way from that of a pseudoplastic fluid; this phenomenon is rather uncommon and is mostly observed in sand-water mixtures. 10In the case where n = 1, the formulation reduces to that of a Bingham fluid.Finally, a bilinear fluid which can be described as a generalization of a Bingham fluid 6 is given by where γ 0 = T 0 /(µ h µ l ).The existence of two different dynamic viscosities µ l and µ h highlights that when the shear rate is high the flow behaves as a Bingham fluid with low viscosity (µ l ), but at low shear rates, the flow behaves as a Newtonian fluid with high viscosity (µ h ). 6he above rheological models have been commonly used for the numerical modelling of cohesive submarine flows. 5,6,11mran et al. 6 have tested the suitability of the three models and observed that each model influenced differently the run-out length of the flow, for example, a bilinear model resulted in a larger run-out distance than a Herschel-Bulkley model.The authors concluded that the bilinear model gives a better representation of the stress-strain relationship, particularly at low shear rates. 6The sediment composition can also play an important role for the choice of the rheological model.The Herschel-Bulkley and especially the Bingham rheology are very common for the modelling of muddy debris flows. 5,7,12n the other hand, for the description of granular flows, the Mohr-Coulomb failure criterion is often used.The rheology of the failing mass (see also Ref. 13) can be expressed as where c is the cohesion of the material (for a purely frictional material c = 0), σ is the stress acting normal to the surface, ι is the pore water pressure, and φ is the Coulomb friction angle.
The term (σ ι) represents the effective normal stress.The flow resistance to motion is a result of the frictional forces dominant between the particles of the mass and the basal layer of the flow.A linear relationship exists between the shear and the effective stress.Three types of shear stresses can be considered the principal mechanisms driving material downslope for submarine landslides: gravity loading, seismicallyinduced stress, and storm wave-induced stresses. 14Gravity driven flow behavior is often incorporated in the landslide modelling.
Earlier studies on the mechanics of submarine landslides indicated that both viscous and plastic terms become important on the simulation of complex flow motion. 14Norem, Locat, and Schieldrop 15 incorporated both visco-and plasticproperties in the rheological regime.The retarding stress could be thus computed (see also Refs. 2 and 14) by where r ι is the pore pressure ratio equal to ι/γh, ι is the pore pressure, and γ = ρg is the specific weight of the sediments with h thickness.The approach in Ref. 15 stresses the importance of incorporating a drag term to account for the forces exerted by the seawater on the frontal and upper body of the flowing mass.The incorporated water drag can reduce significantly the maximum velocity of the landslide.
To account for viscous drag forces, Hutter and Nohguchi 16 made use of a viscous Voellmy-type resistive stress term, proportional to the square of the landslide velocity.More recently, Løvholt et al. 17 have combined a turbulent term with the granular model introduced by Savage and Hutter 13 to incorporate the water drag in their submarine landslide simulations.The turbulent term results from the empirical formula first described by Voellmy for the rheology of snow avalanches. 18The formula can be written as where T is the stress vector, u is the flow velocity, and ξ is a non-dimensional coefficient used to represent the effect of turbulence and/or grain collisions. 19s submarine flows are almost impossible to observe directly, the rheological properties of the flow can be generated by geotechnical measurements in the laboratory.1][22] However, the values of the parameters may not be always easily defined due to the lack of data/facilities, etc.In the absence of adequate input data, a common practice in the numerical modelling of submarine slides is to try to match the depositional lobes resulting from the simulations, with the landslide deposits observed in the field.A drawback in this approach is that the simplifications of the physics in the models cannot capture the complex physics of the landslides.
For example, a constant retarding stress can be often assumed for the whole duration of the flow whereas this is rarely the case in reality.Mixing of the slide material with water increases its pore water pressure and can have a significant reduction in the shear strength of the flow. 5Measurements of the yield strength on the sediments of the slope and on the deposits can thus differ substantially. 23In addition, mechanisms such as hydroplaning and remoulding that enhance the capability of the flow to travel large distances are often not accounted in the modelling.
Lastras et al. 24 have found that for the simulations of the Big'95 debris flow, the best fit solution had a very low yield strength of 800 Pa.The authors have also noted that the internal friction had also a dramatic influence on the run-out length of the flow and best fit values were as low as ∼1 • .They concluded that the absence of hydroplaning in the modelling could explain the very low yield strengths required.Simulations of hydroplaning and non-hydroplaning debris flows in Ref. 7 show that under a non-hydroplaning regime the yield strengths have to be significantly reduced for the flow to reach the same run-out distance.
The sediment composition of the slide material has an effect on the cohesion of the flow.Sediments with high clay content suggest more cohesive landslides and can thus be described by higher yield strengths (such is the case for sliding in the Bear Island Fan). 11On the other hand, mixtures of clay and sand gravel appear to be less cohesive (like the 1929 Grand Banks slide 11 ).Increase in the yield strength of the flow is linked with decrease in the observed run-out distance.De Blasio et al. 5 have observed that an increase in the yield strength from 6 to 14 kPa in the simulations resulted in a decrease in the run-out distance of the flow by a factor of 2. The value of the yield strength can thus be of great importance for the modelling of submarine landslide deposition, and it may vary between different numerical models.
In the present work, we model the episodes of collapse at the Rockall Bank Slide Complex (RBSC), NE Atlantic Ocean under a rheological framework.We note that the terms landslide and slide in this case refer to general slope collapse in the region without implying any particular type of motion.For the landslide simulations, the one-fluid version of the two-dimensional code VolcFlow has been employed. 9olcFlow is a finite-difference Eulerian code that solves the general shallow water equations of mass and momentum conservation.The algorithm can account for a simplified representation of varying rheological regimes.It has been efficiently used to study many cases of subaerial and submarine sliding. 19,25,26 frictional model, based on the Mohr-Coulomb failure criterion, was used to study the flow motion.A simplified Bingham rheology is also implemented.The results of our simulations show that the Bingham rheology represents successfully the landslide deposition for this case, whereas the basal friction critically affects the run-out length of the flow and thus has to be significantly small.To study separate failure episodes, we make use of a model that combines the Bingham rheology with frictional terms for the rheology of the landslide.In all the scenarios, a velocity dependent term is incorporated in the models to consider the effect of hydrodynamic drag in the motion.We make use of a Cartesian coordinate system, where the xand yaxes represent the EW and NS horizontal directions, respectively, and the zaxis represents the vertical direction.The retarding stress T = T x , T y is given by where du/dh is the shear rate, r is the slope curvature, and φ bed is the basal friction angle.As it can be observed, the first three terms of the equation resemble the model introduced by Norem et al. to describe the mobility of subaqueous flows. 15The last term of the equation is making use of the coefficient of turbulence; in VolcFlow, the term is non-dimensional. 19In this paper, we first present an overview of the Rockall Bank and multiphase failure in the area in Sec.I B. The methodology of the work (Sec.II A) and the results of the varying landslide scenarios (Sec.II B) are given in Sec.II.A short discussion on the tsunamigenic consequences of the landslide scenarios is also given in this section (Sec.II C).Finally, the conclusions of this work are presented (Sec.III).

B. Multiphase failure in the Rockall Bank, NE Atlantic
The Rockall Trough is a steep-sided sea-floor depression approximately 1000 km long and 200-300 km wide.It forms the largest of the three basins located northwest of the UK and Ireland in the Atlantic Ocean, the other two being the Hatton Basin and the Hatton Continental Margin 27 (Fig. 1).It is also the closest of the three in proximity to the Irish coast.
The trough is bounded on the east by the Irish and Scottish continental shelves; Rockall Bank, an almost flat-topped, underwater plateau lies to the west (Fig. 1).The bathymetry of the Rockall Trough ranges between 300 and 400 m on its east and west margins, abruptly going down to 3000 m at the north part of the basin, deepening gently to 4000 m toward the south, at the Porcupine Abyssal Plain (Fig. 1).
Two orthogonal wide-angle seismic profiles, recorded during the Rockall And Porcupine Irish Deep Seismic (RAPIDS) project (1988-1990) together with previously attained gravity data revealed that the Rockall Trough formed as a result of a multiphase rift evolution that took place in the Triassic, Late Jurassic, and Early Cretaceous. 27In more detail, the data revealed the existence of continental crust Three major depositional areas are distinguished on the seafloor of the Rockall Trough: the north part is occupied by the Donegal-Barra Fan and the deposits of the Rockall Bank Slide Complex (RBSC), also referred to as the Rockall Bank mass flow (Fig. 1).The areas lie offshore the Irish and Scottish continental shelves.The third one, the contouritic Feni Drift, extends on the west margin of the Rockall Trough and up to the deepest part of the basin in the south 29,30 (Fig. 1).Among the three, the RBSC constitutes the largest region of submarine slope failure in the Irish Atlantic margin, affecting a total area of 18 000 km 2 . 29oberts was the first to describe the submarine failure scarps that extend widely on the eastern slopes of the Rockall Bank. 31Several bathymetric surveys have been carried out in the area ever since.The Geological LOng Range Inclined Asdic (GLORIA) side-scan sonar survey (the Atlantic Irish Regional Survey-AIRS96) took place in 1996 and the Towed Ocean Bottom Instrument (TOBI) survey of a higher 5-10 m spatial resolution in 1998. 32,33Between 2000 and 2001, R. V. Bligh collected multibeam bathymetric data of most of the Irish offshore area as part of the Irish National Seabed Survey (INSS) program.The multibeam bathymetry and sidescan sonar data acquired during the surveys have provided a full coverage of the slope failure features.
Two distinct areas of failure escarpments (upper and lower slope regions) have been identified on the eastern margin of the Rockall Bank. 28Due to its diverse scarp morphology, the upper slope was further subdivided into three distinct parts: the north, the south, and the central upper slope. 28The two regions are separated by an elongate deep erosional moat near the base-of-slope, which traverses the area for 120 km and in parts can be as deep as 150 m. 28,29 The gradients of the slope are steep varying between 5 • and 10 • and in places can reach up to 15 • -20 • . 28The scarps vary on sizes and heights.Despite some isolated retrogressive failure episodes, an overall non-retrogressive failure mechanism has been proposed based on observations of the scarp morphology and geometry. 29ownslope of the scarp complex is the depositional area occupied by sedimentary lobes (Fig. 2).The deposits of the RBSC show a run-out length of 120 km and a depositional width of 150 km 28,29 [Fig.2(c)].The width to length aspect ratio may be characterised as irregular when compared to other submarine landslide deposits found on the western and eastern margins of the Atlantic Ocean. 29,34The occurrence of turbidites in piston cores, collected from the deepest part of the Rockall Trough basin, supports that the flow might have transformed into a more-dilute suspended mixture during motion that had travelled as turbidity current. 35Sedimentary core data suggest a mixture of clayey muds, silt, and sandy layers for the composition of the flow. 35n the past, there has been a high degree of uncertainty regarding the triggering mechanism of the failure and whether sliding occurred as a single failure episode or multiple successive or separate episodes of failure. 28One hypothesis for the triggering mechanism suggests that sedimentation on top of abrupt pre-existing fault scarp topography at the upper part of the slope enhanced by strong bottom current activity removed support from the base-of-slope, leading to slope instability and subsequently failure. 28,29Elliott et al. concluded that the contourite drift dynamics in the basin acted as a substantial long term factor promoting slope instability. 29The sedimentation pattern has caused overloading of the upper part of the slope, whereas the contourite drift systems contributed to the erosion of the lower part of the slope and led to the formation of an elongate moat.Georgiopoulou et al. found a direct relationship between the RBSC seafloor scarps and the basement morphology and suggested that focused fluid flow along basement-bounding faults and/or differential compaction across the scarps contributed to the sliding. 28 manual reconstruction of the preslide morphology estimated that the volume of the missing sediments ranges between 265 km 3 and 765 km 3 (conservative and generous approach, respectively, in Ref. 28).Sedimentary and seismic data, acquired from the depositional lobes of the RBSC during the RV Celtic Explorer cruise CE11011, support the hypothesis that the RBSC morphology is the result of multiple phases of slope failure. 35At least three major phases of collapse can be distinguished based on the data. 28,35In more detail, the seismic profiles along the lateral extents of the lobes have revealed the limits of three depositional horizons separated by long periods of slope stability 35 (Fig. 2).The two first episodes of collapse (here referred as landslide A and B) have more likely initiated from the south upper slope (SUS) region of the margin, whereas the most recent episode (landslide C), the lobes of which are still visible on the seabed of the Rockall Trough, has resulted from collapse on the north upper in combination with failure from the lower slope region. 35s landslide A forms potentially the oldest failure episode and is buried deeply, its age is hard to estimate; it is suggested that it must be a few Ma old. 35An estimation of the average sedimentation rate on the slope scarps suggests that landslide B may have happened at ca. 200 ka. 35Interpretation of radiocarbon dating shows that the most recent event took place 21 ka, which is 3000 yr post Last Glacial Maximum (LGM), when the British Irish Ice Sheet (BIIS) was still at maximum extent and was only just starting to destabilise. 35he volumes of the three phases (landslides A, B, and C) were estimated based on information from the sediment piston cores and the 2D seismic profiles of the landslide deposits at the complex. 35The volume of landslide A was estimated at approximately 200 km 3 , the deposits cover an area of 7500 km 2 [Fig.2(a)]. 35The second phase (landslide B) was smaller in size occupying an area of 4500 km 2 with a material volume of ca.125 km 3 [Fig.2(b)]. 35The data regarding the 3rd failure episode (landslide C) yield a volume of ca.400 km 3 , occupying an area of 6600 km 2 [Fig.2(c)]. 35ollowing these volume estimations, the total volume of the deposits becomes 725 km 3 .This comes in close agreement with the estimated maximum volume of the missing sediments from the slope (765 km 3 ), i.e., the generous approach for the reconstruction of the preslide morphology. 28It is noted that the two older episodes (landslides A and B) involve higher levels of uncertainty than the most recent one (landslide C).The estimations of the depositional extents of landslides A and B are bounded by the limited number of the seismic profiles acquired during the surveys.The event for which the most information exists is the 3rd episode of collapse, the extent of which is visible on the seafloor of the trough.As a consequence, the landslide modelling of the events draws information from the latest episode of collapse to study the landslide deposition.

A. Methodology
For the landslide modelling, we follow closely the methodological approach described in Ref. 36.The bathymetry data were retrieved from the EMODnet Bathymetry portal.The data were transformed to the Cartesian coordinate system, smoothed, and interpolated to form a uniform grid.A Digital Elevation Model (DEM) with dimensions of 200 km × 250 km and spatial resolution ∆x = 320 m, ∆y = 500 m was employed for the landslide modelling.Vol-cFlow computed the run-out length of the landslide.To establish the deposition, the landslide motion was initially modeled with a time step that varied in accordance with the stability criteria in the code.A stopping criterion for a velocity lower than 4 m/s was used, as no significant depositional re-arrangement was observed beyond that point.
For each scenario, a two-dimensional Gaussian function is used to represent the initial body of landslide.The Gaussian shape forms a simplified yet efficient and standardised approach to approximate the landslide body in numerical modelling.8][39] In the numerical simulations of submarine sliding in the RBSC, a landslide of a Gaussian shape was computed with the mass spreading along the xand yhorizontal directions (Fig. 3).The shape of the mass is given as a function of the landslide thickness h and the standard deviation of the landslide, σ x and σ y .
Note that Sarri, Guillas, and Dias 40 performed a sensitivity analysis of the Gaussian shape impact on the tsunami wave created.In that study, there was no significant sensitivity to the spread ratio (i.e., the ratio of the characteristic length over the characteristic width for a constant volume).However, this may be due to the smooth aspect of a Gaussian shape: other irregular (and time-varying) shapes may have a different influence on the resulting wave.However, the investigation of the influence of a high dimensional shape as input is complex and was only recently addressed 41 in the context of the influence of bathymetry onto tsunami wave propagation.
The reconstruction of the most recent episode (landslide C) was performed as more information exists for its depositional lobes. 35Taking under consideration the best-fit values of the rheological parameters for landslide C, the numerical simulations for landslides A and B were performed.This modelling approach follows the hypothesis that the sediment composition in the area is fairly similar.As introduced in Sec.I B, the slope failure that has resulted in the depositional lobes of the last episode possibly initiated from the North Upper And Lower Slope (NUS-LS) regions.The scarp heights on the slope region are measured to range between 200 and 300 m at the north upper and south upper slopes.The lower slope is dominated by smaller scarp heights that range between 60 and 120 m.
These scarp heights were used for the initial thickness of the landslide.A landslide with maximum thickness h max NUS = 300 m and a volume of V NUS = 264 km 3 was assumed for the north upper slope.For the lower slope region, a landslide of h max LS = 120 m and volume V NUS = 136 km 3 was considered (Fig. 3).We took under consideration the maximum scarp heights of the landslides first to account for a worst case scenario and second to avoid excessively large spatial extents of the landslide on the slope.Finally, the total volume of the event was V 3 = 400 km 3 matching the estimated volume of landslide C. For the older episodes, we consider two landslides initiating from the south upper slope region (landslides A and B).Landslide A corresponds to the oldest recorded failure episode in the region with an estimated volume of ca.200 km 3 . 35Another collapse initiating from the same slope area with a volume of approximately 125 km 3 (landslide B) happened at a much later stage. 35

B. Landslide scenarios
Table I presents landslide scenarios examined with the three rheological models.Seven landslide scenarios are presented here, of which scenarios 1-5 simulate landslide C and scenarios 6 and 7 simulate landslides A and B, respectively (Table I).Scenarios 1 and 2 make use of a purely frictional rheology similar to Eq. ( 4), scenarios 3 and 4 represent a Bingham rheology with frictional terms [described by Eq. ( 7)], and scenario 5 studies the effects of the Bingham rheology on the landslide motion [Eq.( 1)].The main rheological parameters incorporated in the numerical simulations are the yield strength The one-fluid version of VolcFlow assumes that no mixing is taking place between the water and the landslide.To partially account for mixing, a reduced density is used for the landslide (ρ = ρ ls ρ , where ρ ls denotes the landslide density and ρ denotes the water density, respectively).In numerical modelling of submerged motion, such an approach has been used for many cases. 20,26,42Here, we consider a fixed density of ρ = 1, 200 kg m 3 .
Landslide scenarios 1 and 2 model a frictional flow motion with the incorporation of a velocity dependent term [Table I and Figs.4(a) and 4(b)].It is observed that the basal friction has to be very low to produce a run-out length [indicated with a black dashed line in Fig. 4(a)] that matches with the one of the field deposits (orange dashed line).In more detail, a flow that has φ bed > 0.5 • results in deposits that do not move significantly but rather accumulate close to the slope region.They exhibit a more lateral than longitudinal spreading [Fig.4(a)].The basal friction has to be significantly reduced down to φ bed = 0.2 • to produce a run-out length that matches the observations [scenario 2, Fig. 4(b)].However, in this case, the landslide deposits exceed the lateral limits of the depositional lobe and also exhibit large spreading close to the slope region [Fig.4(b)].
It is noted that the dimensional coefficient ξ in both scenarios is equal to 0.045 (Table I).The term has a direct impact on the maximum velocity of the landslide which rises up to 60 m/s.Maximum velocities of such values are not unlikely for energetic landslides. 43However, as the validity of large maximum velocities is under debate (see also discussion at Ref. 36), the value of the drag is not further decreased.An increase in ξ results in smaller flow velocities.The deposits of scenario 3 are computed with ξ = 0.055, for which the maximum velocity of the landslide becomes approximately 50 m/s.
The value of the maximum velocity falls in agreement with values of similar events. 36The simulated deposits of the scenario reach the run-out length of the deposits observed on the field [Fig.4(c)].Albeit, the computed deposits extend significantly beyond the lateral limits of the observed lobes in the NE-SW direction.This is possibly linked to the low yield strength of the sliding material (T 0 = 550 Pa).Indeed in scenario 4, we increase the yield strength to T 0 = 750 Pa which results in more cohesive/constrained flow deposits [Fig.4(d)].In this case, however, ξ has to be decreased to get an adequate run-out length (ξ = 0.045, Table I).The decrease in the drag term allows for an increase in the value of T 0 .In this scenario, the simulated deposits match better the observations [Fig.4(d)].
In scenario 4, the dynamic viscosity of the flow is also tuned to 300 Pa s, similar to the dynamic viscosity of flows with large extent as discussed in Ref. 20 (Table I).This change in the values of viscosity does not have any significant effect on the run-out length of the flow (scenario 3 compared to 4 [Figs.4(c viscosity up to two orders of magnitude in the simulations do not have significant effect on the results. 20In Ref. 36, it was also observed that the best values of µ exhibit a uniform distribution. Landslide scenario 5 simulates the landslide as a flow that exhibits behavior similar to that of a Bingham fluid but with the incorporation of a drag term.In that case, the yield strength of the landslide can be increased up to 1500 Pa for a good representation of the observations.The values of the other input parameters are kept constant (scenario 5, Table I).The computed maximum velocity of the landslide exhibits a similar trend to that of scenario 4. The simulated run-out length of the deposits falls in good agreement with the field deposits [Fig.4(e)].The lateral extent of the simulated deposits of scenario 5 is smaller than the extent of scenarios 3 and 4 (slightly less than scenario 4).The flow appears to be more constrained and thus in better agreement with the field observations.
From the scenarios examined, scenarios 4 and 5 form the best-fit cases for the Bingham and the Bingham with plastic terms regime (Table I).The "extended" Bingham model is possibly the best approach to simulate landslide C with VolcFlow.From a numerical point of view, increasing the basal friction of the landslide does not have any beneficial effect in the results.As the yield strength has to be reduced for the landslide to reach the desired run-out length, the sliding material spreads laterally and gives wider deposits than observed in the field (as exhibited in scenarios 3-5).However, the side scarps observed on the lower slope region constitute strong indicators of erosion and therefore dragging during the flow motion.This phenomenon renders the presence of basal friction nonnegligible (with values that in reality could have been higher).For this reason, we choose the values of scenario 4 to simulate flow motion from the SUS region.
The simulated deposits of landslide A reach a run-out length similar to that of landslide C [Fig.4(f)].These results come in general agreement with the estimations of the extent of landslide A in the field. 35However, the simulated lateral extent and run-out length of landslide A (black dashed line) are somewhat smaller than the estimated extent of the field deposits (yellow) [Fig.4(f)].It could be possible that the landslide in reality was more dilute than landslide C and the deposits have spread more.In the seismic profiles of the landslide deposits, the thickness of landslide A is found much smaller (∼30 m) than the thickness of landslide C. 35 Simulations of landslide A with a lower yield strength (T 0 = 550 Pa) and similar rheological properties, as in scenario 3, would result in deposits extending slightly more laterally.
In scenario 7, the smaller volume of the landslide (125 km 3 ) acts as a constraining factor for the run-out length [Fig.4(g)].In more detail, the toe of the simulated flow comes in rest before it reaches the run-out length of landslide B (purple line) [Fig.4(g)].In this case, the simulated deposits spread more laterally than the field deposits.It is noted, however, that due to the limited number of available seismic profiles for the area, the estimation of the landslide deposits for landslides A and B forms an approximate interpretation of the data and in reality the deposits could have had a different extent or shape.Additional research in the area would be necessary in order to obtain more insight and advance the numerical modelling of those events.
The contours of the landslide thickness after deposition are also presented for scenarios 4-7 (Fig. 5).The thickness of the deposits for scenarios 4 and 5 appears to be fairly similar with a higher concentration of sediments closer to the landslide toe.Scenario 4 results in deposits of slightly higher thickness compared to scenario 5 [Figs.5(a I).The colour-scale indicates the thickness (m); the numbers on the deposits correspond to the thickness measurements (15, 30, 45, and 60 m).match with the estimations of the real landslide thickness from the seismic profiles, which is for both cases estimated at approximately 30 m. 35 Overall, it is noted that for good agreement of the simulated run-out lengths with the field observations the values of the input parameters have to be significantly low.The low or negligible basal friction that is required for the flow to reach this run-out length may be attributed to complex processes that the code cannot capture.The same applies for the yield strength of the material.Such are the changes in the rheological composition of the mixture during motion, (especially when water is mixed with the flow) which are not uncommon in submarine sliding.Bradshaw, Tappin, and Rugg, 23 for example, have found values of T 0 = 30, 38, and 45 kPa on slope sediments, whereas they measured T 0 = 1, 10, and 20 kPa, respectively, on the deposits of the debris flow. 23ther physical processes are the incorporation of water in the flow, which is only partially considered in this study.Processes that involve more uncertainty are the addition of material through erosion (which could result in volume changes during the motion) and the possible existence of an element of regression that is not visible in the current morphology.This could have an effect on the slope failure mode and possibly tsunami generation.It is also not known yet whether these events occurred independently or if they constitute composites of smaller events (i.e., comprising of smaller individual slides).
The significantly low yield strength of the modeled flow in Ref. 24 was linked with additional flow mechanisms such as remoulding and hydroplaning that can enhance the flow mobility.Hydroplaning is strongly linked with the cohesiveness and the sediment composition of the landslide. 20Less cohesive flows with a sand-rich sediment composition can reach long run-outs without hydroplaning, whereas clay-rich flows tend to be more cohesive, and thus they are more likely to hydroplane. 11Hence, high yield strengths may indirectly lead to greater run-out distances.However, the possibility of hydroplaning taking effect in the RBSC is not known with certainty.The erosive nature of the collapse may render hydroplaning an unlikely process for the region.

C. A glimpse on consequences
The nature of the failure mechanism and the large volume of the failed sediments raise questions for the tsunamigenic potential of the failures, as they form factors that can promote tsunami generation. 44As shown in Ref. 36, submarine failure in the RBSC has a tsunamigenic potential that could affect the northwest Irish shoreline.As the age of most recent failure episode (landslide C) falls between 23 and 19 ka BP, 35 any landslide tsunamis generated in this periods would probably not have affected Ireland since the western shelf and coast were still covered by ice. 45In more detail, the reconstruction of the different retreat stages of the BIIS and the North Sea ice cover shows that the entire country was covered by the BIIS 23 ka ago. 45At 19 ka, some parts of the northwest coastline and shelf might have been revealed; nonetheless a significant part of the west coast was still under ice. 45There is a high likelihood that any propagating tsunamis generated in this period would have reached an ice-covered shelf.This current hypothesis excludes the possibility of locating any tsunami deposits on the northwest Irish shoreline and therefore the validation of current and future modelling results.
It is noted that the uncertainty surrounding landslide C in the RBSC is large.Even if the above hypothesis is true, assessing the risk of submarine failure in the RBSC taking Phys.Fluids 30, 030705 (2018) FIG. 6. [(a)-(c)] The free surface elevation is captured by the offshore wave gauges 1-3 (at water depths of 2632, 2838, and 1383 m, respectively).The results of different scenarios are represented with different colours.In more detail, scenario 3 (blue), scenario 4 (orange), scenario 5 (yellow), scenario 6 (purple), scenario 7 (green) (see also Table I).(d) The computational domain and the locations of gauges 1-3.
as an example the current topography of the region can be of great importance to understand the implications that a similar event could have if it occurred now.A comprehensive tsunami hazard assessment for the region is currently left for future research.In this study, we show the results of some simulations with the Nonlinear Shallow Water Equation (NSWE) solver VOLNA. 46VOLNA gives a simplified overview of tsunami generation and propagation in the region without accounting for the complex characteristics of landslide tsunami propagation.The numerical approach of the computations (e.g., computational grid, coupling with VolcFlow) is discussed in Ref. 36.The GEBCO 08 GRID terrain model was used for the land elevation.
We present here three offshore wave gauges that measure the free surface elevation over time at water depths of 2632 m (G01), 2838 m (G02), and 1383 m (G03) [Fig.6  also given in Table II.The results show that the largest wave amplitudes are produced by landslide C.Although the initial free surface elevation between landslide A and B does not differ significantly (Table II), the discrepancy between the two scenarios increases as the tsunamis propagate [Figs.6(a)-6(c), Table II].
Because of its lower kinetic energy, scenario 3 results in smaller tsunami amplitudes than scenarios 4 and 5 [Table II, Figs.6(a)-6(c)].In fact, scenarios 4 and 5 result in tsunami amplification of similar magnitude, which indicates that the rheological regime in the examined case does not affect significantly tsunami amplification.This comes in agreement with the results of rheological parametrisation studies. 21,38The free surface elevation generated by scenarios 4 and 5 gives similar records on all the gauges with scenario 4 generating tsunamis of slightly larger amplitude than scenario 5 (Table II).A more cohesive landslide could thus result in slightly larger tsunami amplitudes.However, general conclusions about the relationship between the parameters cannot be drawn by this small sample.Furthermore, this difference is less than 5% of the tsunami amplitude of the largest scenario.Given the amount of the uncertainty that underpins the event, these discrepancies are insignificant.
The first waves that arrive in gauge 01 are the waves from scenarios 3, 4, and 5, probably due to the location of the gauge closer to the NUS-LS region [Fig.6(d

III. CONCLUSIONS
Overall, it is exhibited that the varying geomorphological and rheological landslide input parameters may affect significantly the depositional process of the landslide and the tsunami wave generation.Careful consideration has to be given when addressing the issue of the landslide prefailure and post-failure conditions.Different rheological models might result in diverse landslide deposits but similar tsunami wave characteristics.On the other hand, the same rheological parameters may result in diverse wave characteristics when the geomorphological characteristics of the landslide differ substantially.
To counterbalance the lack of the water drag in the onefluid version of VolcFLow, a drag term was introduced.In these scenarios, the peak velocities of the flow were in good agreement with the current knowledge for similar events; however, the resistance of the flow to motion was enhanced.As a result, other parameters, like the basal friction and the yield strength of the sliding material, have to be significantly reduced in order to compute the observed run-out length of the flow.Lower values of these parameters contribute to the lateral spreading of the deposits.Additional mechanisms that may enhance the flow mobility and reduction of the yield strength during the motion were not modeled by the code; these phenomena are also difficult to record in nature due to the scarce direct observations of submarine sliding.
The best-fit solutions for landslide C were obtained with the use of the Bingham rheology including a velocity dependent term (scenario 5, Table I), or the use of a rheology that incorporates the Bingham and frictional flow properties in the motion, including again a drag term (scenario 4, Table I).Comparisons of both regimes show that the Bingham rheology results in better constrained landslide deposits; however, incorporation of basal friction may be closer to the true nature of the failure.The tsunami simulations of the last episode of multiphase collapse in the RBSC with both regimes yield very similar tsunami amplitudes, albeit slightly higher in the case of a flow without friction [Figs.6(a)-6(c)].Landslide simulations with a frictional rheological model were also examined but were not considered advantageous for this case.
Following the best-case scenario numerical modelling of the older landslide episodes (A and B) in the south upper slope region was performed.The modeled deposits do not match as well with the observations.This can be attributed either to submarine sliding with different rheological properties (e.g., a more dilute flow) or to the current lack of knowledge regarding those events due to data limitations.Additional field and numerical research will be beneficial for a more comprehensive analysis of these episodes.In addition, a more thorough statistical approach for calibration could be used to better tune the input parameters.However, the dimension of the outputs is large since ideally both the extent curve and the thickness surface have to be compared to observations.A step in this direction has been taken recently to reduce dimension without sacrificing accuracy in the calibration. 47The ultimate challenge would be to model, and thus calibrate, the timevarying change in the rheological composition of the mixture (and in volume assuming deposits) during motion, but these are currently intractable questions.
The landslide scenarios examined show that slope failure in the region has given rise to tsunamis of varying magnitudes.After tsunami generation, a leading tsunami crest is formed that propagates toward the Irish coastline followed by a tsunami trough.The seabed topography of the Rockall Trough strongly influences the tsunami wave propagation.This can be especially observed when the waves reach the east margin of the Trough.Reflection and refraction occur above the continental shelf as the waves pass from one medium to another.The first area to be inundated in the Irish coast is the Mullet Peninsula, in Co. Mayo.The time arrival of the first tsunami wave peak in the peninsula may be less than 50-55 min after landslide generation.At the middle part of the peninsula, tsunami overtopping is observed.The waves enter Blacksod Bay in two ways: overtopping of the land in the middle of the peninsula and propagation from the entrance of the bay.
The tsunami amplitudes and direction depend on the scenario examined (volume, maximum thickness of the landslide, region of collapse).Focusing on the three major episodes of multiphase collapse attributed to the RBSC, the most voluminous event (landslide C) has the largest tsunamigenic potential.Gauges close to the generation region measure a maximum tsunami amplification of approximately 40 m that decays significantly during propagation.Landslide-induced tsunamis from the SUS region are smaller in amplitude.The second episode of collapse (landslide B) holds the smaller tsunamigenic potential of the three.
The results of this study constitute an attempt to study numerically the different episodes of collapse in the RBSC under a rheological framework and offer a view on subsequent tsunami generation.Future research could be devoted to a more sophisticated modelling of the phenomena.This could be approached by employing a fully coupled numerical code.In addition, to adequately assess the tsunami propagation and predicted run-up heights on the coast, the use of more advanced dispersive solvers can improve the modelling in some settings.Finally, a slope stability assessment at the source and a different meshing strategy in the modelling with mesh refinement on the coastal areas (where impact would be better estimated as a result) would be of great value for a comprehensive hazard assessment in the area.

FIG. 1 .
FIG. 1. Bathymetric map of the Rockall Trough with the main morphological features highlighted.The bathymetric data are based on GEBCO and the contours are in meters.The white box indicates the area of study.The stars indicate roughly the positions of the slope regions (north upper slope-red, south upper slope-yellow, lower slope-black).For a detailed description of the regions, we refer to Ref. 28.

FIG. 2 .
FIG. 2. Shaded relief bathymetric map of the RBSC showing the distribution of landslides A, B, and C as derived from the interpretation of the sedimentary core data, the seismic profiles, and depositional lobes on the seafloor.Bathymetric data from the Irish National Seabed Survey.From Ref. 35.

FIG. 3 .
FIG. 3. (a) The slope region in 2D with the integration of landslide C. The red star shows the landslide centre on the north upper slope region, and the black star shows the landslide centre on the lower slope region.(b) The 3D plot of the slope region.The colour-scale shows the depth in meters.
FIG. 4. Simulations of the landslide deposits deriving from scenarios 1-7 (a-g respectively) in Table I.The black dashed lines indicate the numerical depositional extent.The red dashed line in (a)-(e) exhibits the extent of the deposits for landslide C. The yellow and purple lines in (f) and (g) show the estimated depositional extent of landslides A and B, respectively.
FIG.5.The thickness of the numerical deposits of scenarios 4-7 (a-d respectively) (TableI).The colour-scale indicates the thickness (m); the numbers on the deposits correspond to the thickness measurements(15, 30, 45, and 60 m).
(a)].The gauges record the tsunami amplitudes close to the source region and during tsunami propagation over the Rockall Trough; the tsunami simulations for scenarios 3-7 are presented [Figs.6(a)-6(c)].The initial water surface elevation is FIG. 7. Snapshots of the tsunami propagation: (a) over the Rockall Trough at t = 45 min.The black box shows the area in b.Black arrows indicate reflection waves when the tsunami passes over the edge of the Rockall Trough (b) close to the Mullet Peninsula in Co. Mayo at 65 min.WSE refers to the Water Surface Elevation in meters.

Figure 7 (
b) shows the tsunami inundation with a smaller colour scale to increase the intensity of the inundation; it can thus be observed that the waves overtop the middle part of Mullet Peninsula.At t = 65 min, the tsunamis enter the bay in two ways: the one is by overtopping of the peninsula and the other one is by the opening between the Mullet Peninsula and Achill Island [Fig.7(b)].Merging of the two tsunami wave peaks and tsunami oscillations are likely to occur inside the bay at a later stage.The Inishkea Islands seem to act as barriers for the propagating waves, reducing the tsunami amplification behind them [Fig.7(b)].Reflection of the waves is also observed on the north part of the peninsula, the Inishkea Islands and Achill Island.After two hours of propagation, the waves have reached the Aran Islands in the south.

TABLE I
T 0 , the basal friction angle, φ bed , the dynamic viscosity, µ, and the coefficient of turbulence, ξ, used to represent the hydrodynamic drag.The selection of the appropriate values of the rheological parameters is based on the literature and the results of the Bayesian calibration in Ref.36.

TABLE II .
Table II illustrates the initial water surface elevation (WSE) and the recorded surface elevation in gauges 01-03 [Figs.6(a)-6(c)] in meters.