Modeling the hydraulic excavation of cohesive soil by a moving vertical jet

Abstract In dredging, high pressure water jets are commonly applied to assist the mobilization of soil. This work considers the excavation of cohesive soil. The key objective is to predict the development of the cavity in the soil as a function of the undrained shear strength, translation velocity and hydrodynamic pressure of a single nozzle. A generic computational fluid dynamics (CFD) model has been developed that captures both the jet flow and the soil failure in a single framework. The results are compared with data from a previous experimental study. The CFD model predicts the cavity dimensions with reasonable accuracy. In addition the model provides detailed data to study the cyclic nature of the soil failure process. The CFD model is promising and can be applied for more complex nozzle configurations to assist the design process of dragheads and improve production estimates.

This study focusses on the hydraulic excavation by a translating high pressure jet for the excavation of cohesive soil, where the dynamic pressure of the jet, p j , is typically much higher (≫5) than the undrained shear strength su. The soil failure process under these conditions is very complex and the erosion type is classified as mass erosion (Winterwerp et al., 2012). The soil fails undrained and lumps of soil are torn from the bed by the hydrodynamic stresses when the jet traverses the bed (Fig. 2).
Mass erosion is much harder to predict compared to floc and surface erosion that correspond to milder flow conditions in the marine environment (Winterwerp et al., 2012). The case of a fixed submerged jet on cohesive soil has received a lot of attention in literature e.g. (Mazurek et al., 2001;Mazurek and Hossain, 2007;Hou et al., 2016;Dong et al., 2019) and standard jet erosion tests have been developed to assess the erodibility of cohesive soils (Hanson, 1991;Wardinski et al., 2018). Translating jets on non-cohesive soils are also studied frequently (Perng and Capart, 2008;Yeh et al., 2009;Weegenaar et al., 2015). The combination of a translating jet above cohesive soil has received less attention. Zhang et al. (2017) explored the high discharge regime with a relatively low jet pressure in a laboratory experiment (1.2 times su) and a sea trial (1.7 times su). Rockwell (1981) is, to the best knowledge of the authors, the first publicly available report considering the high-pressure regime. He conducted a series of laboratory experiments and explored scaling relations for the cavity depth based on dimensional analysis. Machin and Allan (2011) demonstrated that it is important to distinguish two contributions in the cavity formation: a quasi-instantaneous failure component related to the bearing capacity of the soil and an erosion part. The later contribution strongly depends on the exposure time and hence the traverse velocity of the jet. They also argued that the erodibility depends on the permeability and the plasticity of the soil and concluded this relation requires further study. In their experiments the role of erosion on the cavity formation was considered negligible for a jet traverse velocity higher than 0.1-0.5 m/s. For lower speeds significantly deeper cavities where formed due to progressing erosion at the cavity surface. The model of Machin and Allan (2011) for the maximum cavity depth requires an estimate for the bearing capacity factor (N c ). A constant value corresponding to shallow cone penetration was suggested. However, an appropriate choice for N c was regarded by the authors as an open question. Nobel (2013) conducted a large set of experiments. Different types of failure were identified depending on the relative jet pressure p j /su and the traverse velocity v t of the nozzle. More details about these observations and regime definitions that are relevant for the present study will be provided in Section 1.1. Nobel (2013) argues that both the static pressure and shear stress distribution on the cavity wall are important. The magnitude of each force depends on the location in the cavity. The experiments of Zhang et al. (2016) confirmed the scaling of the cavity depth with the nozzle diameter as proposed by Rockwell (1981). In contrast to the quasi-instantaneous penetration concept of Machin and Allan (2011) and the observations of Nobel (2013), Zhang et al. (2016) argue that erosion at (Yan et al., 2020) the cavity surface controls the migration of the cavity surface and the role of static pressure can be ignored. This assumption would strongly simplify the modeling as it would allow the usage of existing erosion relations for cohesive soil determined for unidirectional flows along the trajectory of the jet. Unfortunately, information about the static pressure at the cavity surface is missing in literature so support for this argument is lacking and requires further study. More recently, experiments of Yin et al. (2019) on extremely high pressure jets up to 30 MPa confirm the general trend reported in previous works of a decreasing cavity depth with increasing traverse velocity and increasing cavity depth with increasing jet pressure. The above literature survey demonstrates that many questions remain open regarding both the qualitative understanding of the failure mechanism and the prediction of the penetration depth and width of the cavity.
Existing analytical models strongly depend on assumptions on the driving forces and failure mechanisms. Computational Fluid Dynamics (CFD) could be helpful as it is a more generic approach and it could reveal diagnostic information that is difficult or even impossible to retrieve experimentally such as the static pressure on the cavity walls and shear plane formation in front of the jet. CFD models have been developed for impinging jet flow and scour related problems Yan et al. (2020). The challenges for a CFD model for this problem are to track or capture the soil-fluid interface (Charin et al., 2017) and to model the undrained cohesive failure in the soil region simultaneously. In this work we follow the methodology proposed by (Lalli et al., 2005) to combine the fluid and solid region in a single computational framework by applying suitable constitutive models for the stress tensors in both regions. The advantage is that mesh deformations are not needed and external coupling between soil mechanical and fluid mechanical solvers is avoided. Recently Wang and Song (2019) have also followed this modeling strategy for static jet excavation in cohesive soil. The key objective is to predict the cavity shape as a function of the jet traverse velocity, jet pressure and the undrained shear strength of the cohesive soil.

Failure modes during jetting
In this work the failure mode definitions of Nobel (2013) are adopted, which will briefly given here for the sake of readability of the present work: • Penetrating jet ( Fig. 4): (1) Penetrating jet happens when p j /su > 12 and v t ≤ u f,h , where u f,h is the soil horizontal propagation velocity, which is defined in Fig. 3. (2) The soil cavity is narrow and deep, with a cavity width of 1-1.5 times the jet diameter. (1) Deflecting jet happens when 7.3 < p j /su < 12 or v t > u f,v , where u f,v is the soil vertical propagation velocity, which is defined in Fig. 3. (2) The soil cavity is shallower compared with the penetrating jet and smaller than 2.5 times the nozzle diameter.
(3) A soil wall texture with curving nerves appears, representing the deflection zone. (4) Limited amount of dislodged soil lumps can be found after the experiments. • Dispersing jet: (1) Dispersing jet happens when 5.4 < p j /su < 7.3 or v t ≫u f,v .
(2) The depth of the soil cavity is shallow while the width of the soil cavity is wider, which can be equal to 5 times the nozzle diameter.
(3) The soil wall structure is irregular. (4) The dislodged soil can be found after the experiments.    (Nobel, 2013).
• Hydro-fracturing: Hydro-fracturing happens when v t < 0.15 m/ s. The soil cavity dimensions are irregular. The soil fails at its weakest position.
The penetrating jet and deflecting jet are the most relevant failure modes. The penetrating jet and deflecting jet are the main objectives to be investigated in this study. The penetrating jet has two jet zones, which are the non-deflection zone and the deflection zone, see Fig. 4. In the non-deflection zone of the penetrating jet the pressure builds up on the top of the soil element until the jet stagnation pressure exceeds the bearing capacity of the soil. Thereafter the soil is pushed downwards by the stagnation pressure of the jet. The smaller the jet ratio (p j / su), the longer the duration of this vertical process. When the traverse velocity of the nozzle is relatively low, the jet flow can follow a vertical line in the non-deflection zone as shown in Fig. 4. As the traverse velocity of the jet is increasing, the jet flow in the non-deflection region will have an inclination angle as shown in Fig. 4. The higher the jet traverse velocity, the larger the inclination angle will be. In the deflection zone of the penetrating jet, the jet flow will be deflected due to the traverse velocity of the jet and the resistance of the soil.
As the jet traverse velocity increases even more, the non-deflection zone will disappear, which is shown in Fig. 5. This failure mode is called deflecting jet. The failure mechanism of the deflecting jet is comparable with that in the deflection zone of the penetrating jet.

Material and methods
Inspired by the works of Nguyen et al. (2014), Weij et al. (2016), Goeree (2018) and Wang and Song (2019) a model for the solid-fluid mixture is combined with a plastic description of soil to capture the jet flow and the cavity formation in a single computational framework. The drift-flux model describes the solid-fluid mixture and the cohesive soil is modelled as a Bingham plastic with a fixed yield stress to represent the undrained shear strength. A threshold value of the solid phase fraction distinguishes solid from fluid regions in order to select locally appropriate constitutive laws for the mixture flow.

Drift-flux model
The drift-flux approach consists of a mass and momentum balance of the solid-fluid mixture combined with separate transport equations for the solid phase fractions.
The mass balance of the mixture results in the continuity equation: where ρ m is the mixture density, u m is the mixture velocity. The mixture momentum equation is: where p m is the mixture pressure, T m is the mixture shear stress tensor. g is the gravitational acceleration vector. α k is the volume fraction of phase k. This term is closed by the transport equation: where Γ t is the turbulent diffusion coefficient. u km , which is the relative velocity between mixture and phase k, is closed by the following equation:

Fig. 4.
Definition of jet zones of penetrating jet. Zone I is the non-deflection zone, zone II is the deflection zone. When the jet traverse velocity is larger, jet flow follows a inclination line (Nobel, 2013).
where c k is the mass fraction of phase k, u kr is the relative velocity between the carrier fluid and the solids which is obtained by the hindered settling relations of Richardson and Zaki (1954). The mixture pressure p m is governed by a Poisson equation that can be obtained by combining the continuity and momentum Equations (1) and (2) respectively, for details see Goeree (2018).

Cohesive soil model -Bingham plastic
To implement the cohesive soil into the drift-flux model, the non-Newtonian fluid model Bingham plastic (Bingham, 1922;Lalli et al., 2005;Goeree et al., 2016) is applied. Bingham plastics behave as a solid at low shear stress and as a fluid at high shear stress which is represented by the constitutive relation: where τ y is the yield stress, η is the plastic viscosity, γ is the shear rate in 3D equal to the second invariant of the rate of deformation tensor ∇u m . From Equation (5), the viscosity of the Bingham plastic is calculated as: where ε is an arbitrarily small value, which is used to prevent singularity in the limit of vanishing shear rate. The yield stress τ y is only applied in the soil region, which is determined by the volume concentration. The yield stress τ y follows: where α s is the volume concentration of the sediment phase, α t is the threshold volume concentration value.

Mixture stress tensor
The RANS model buoyant − k − ε is used in this study to compute the eddy viscosity. In buoyant − k − ε model (Henkes et al., 1991), the modulation of turbulent kinetic energy caused by density stratification effects is incorporated. For the sake of compactness, the equations for the turbulent kinetic energy k and dissipation ε are omitted here. The eddy viscosity is computed as μ e = ρ m C μ k 2 /ε, where C μ = 0.09 is a universal constant. The mixture stress tensor, required in the momentum balance Equation (2) is finally computed as: where I is the identity tensor. The mixture viscosity μ m equals the regularized Bingham viscosity, Equation (6) in the solid region and equals the eddy viscosity in the fluid regions (Boussinesq approximation of turbulent stresses). It is assumed that after dislodgement the soil elements disrupt quickly into fine solid particles (clay) that instantaneously follow the turbulent motions of the fluid carrier phase. For this reason it is assumed that the turbulent diffusivity Γ, required in Equation (3), equals the turbulent viscosity ν e = μ e /ρ m , which corresponds to a turbulent Prandtl-Schmidt number σ = 1.

Numerical implementation
The finite volume method has been selected to solve the governing equations described in the previous sections. All flow variables are defined at the cell centres of the mesh (collocated mesh). The momentum, transport and Pressure-Poisson equations are solved implicitly in a segregated approach inspired by the PISO algorithm of Issa (1986). The open source platform OpenFOAM v1812 (open field operation and manipulation) is used for CFD numerical simulation in this study (OpenFOAM, website). OpenFOAM is an open source object-oriented library for numerical simulations in continuum mechanics based on the finite volume method. In this study the drift-flux solver (drift-FluxFoam) is combined with two dynamic mesh algorithms arbitrairy mesh refinement (A/R) and arbitrairy mesh interface (AMI) to model the translation of the jet.
As for the numerical schemes used in this paper, a central scheme has been applied to the diffusion terms of the momentum equations. The Van Leer scheme has been implemented for the concentration convection terms. The local min scheme, which takes the minimum value of the owner and its neighbour cells, has been implemented to the turbulent eddy viscosity to prevent numerical diffusion at the soil water interface. A threshold concentration for the turbulence model is also included in the model. When the concentration in the cell is larger than the threshold value, the turbulent eddy viscosity will be switched off. By doing this, the turbulence eddy viscosity can be switched off inside the soil.
To model the moving jet, dynamic mesh methods are used in this study. In this study, a combination of AMI (arbitrary mesh interface) and A/R (cell layer addition removal) is used to model the moving jet, see Fig. 6.

Mesh techniques
The steady mesh region separates the A/R dynamic mesh region by the AMI interface. The moving jet is located in the A/R dynamic mesh region. When the jet is moving from left to right, new mesh layers will be added to the left blue region, old mesh layers will be removed from the right yellow region. The jet has a diameter of 30 mm and a SOD (stand off distance: the distance between jet outlet and soil surface) of 20 mm which are representative values for dredging engineering applications (Nobel, 2013). A Cartesian computational grid is applied with a mesh size of 0.005 m.

Initial and boundary conditions
The initial condition is presented in Fig. 7. The height and width of the soil region are 0.28 m and 0.15 m, which ensures the soil cavity is smaller than the computational domain. The side view and front view of the computational domain with boundary conditions setup is shown in Fig. 8. The cohesive soil is placed in between the two side walls. Outlet boundary conditions are implemented in front, end and top of the computational domain. No slip boundary condition is applied to the walls. Dirichlet boundary condition is implemented to the jet inlet. Wall functions (Wilcox et al., 1998) are applied to the turbulence dissipation rate ε and turbulence kinetic energy k on the walls. Table 1 shows the detailed boundary conditions setup of the computational domain. It has to be mentioned that no boundary condition is applied to the interface between water and soil. In this numerical model, the concentration at the soil water interface is not changing sharply but gradually.

Input parameters
A key parameter for the hydraulic excavation process is the ratio between the jet pressure p j and the undrained shear strength of the clay (su) (see Section 1.2). This study focuses on the regime p j ≫su. The second key parameter is the transverse velocity v t . The input parameter matrix covers different combinations of p j /su, traverse velocities v t and nozzle exit velocities v j , see Table 2. The nozzle diameter and stand off distance are kept constant. These settings correspond to the operational regime of nozzles in the offshore and dredging field. The only difference is that the effect of a single nozzle is studied here, while in practice multiple nozzles operate simultaneously and are often combined with mechanical tools. Comparable settings were considered in the experimental study of Nobel (2013). The values for the solid concentration α s and soil mixture density ρ m are adopted from the latter experiments.

Results and discussion
In this section, the results calculated by the CFD model will be validated with the experiment of Nobel (2013). Some analysis of the soil failure mechanism based on the CFD model are also carried out. An adjustable time step is applied based on the Courant number stability criterium. At the same time, mesh and residual convergence analyses are carried out before conducting the validation. Fig. 9 shows the jet cavity during jetting calculated by the CFD model. When the jet traverses the soil sample the cavity formation is clearly retrieved by the model. It can be concluded that the model qualitatively reproduces the general picture of the cavity formation process for high-pressure jets in cohesive soil as sketched in Fig. 2. The jet flow quickly starts to disperse after the deflection, which reduces the velocity considerably.   Note: u is the velocity. p rgh = p − ρgh, which is the pressure without hydro-static pressure. ρ is the local density and h is the water depth. k is the turbulence kinetic energy, ε is the turbulence dissipation rate and α s is the solids volume concentration.

Validation with experiment
3.1.1. Failure modes during jetting Two important jet regimes are shown in Fig. 10. It is clearly seen in Fig. 10 (a), that in Sim 1 two zones are present: the non-deflection zone and the deflection zone as defined in Section 1.1. Therefore it can be concluded that the model can represent the penetrating jet regime. Fig. 10 (b), shows that the jet flow immediately bends backwards indication that only a deflection zone is present. This means that the model can retrieve the regime change between the penetrating jet and deflecting jet regimes.

Cavity depth and width
The dimensionless soil cavity depth is calculated as Z c /D n , where Z c v j is the nozzle exit velocity, p j /su jet pressure divided by the undrained shear strength, v t is the jet traverse velocity, D n is the nozzle diameter, SOD is the stand off distance, α s is the solid concentration, ρ s is the solids density.  is the soil cavity depth and D n is the nozzle diameter. To validate the CFD model, CFD Sim 2, which has a jet ratio p j /su = 19.5 with different jet traverse velocities v t , has been carried out. The results are validated with the experiment of Nobel (2013), see Fig. 11. The buoyant − k − ε RANS model has been applied and compared with the experimental results.
The uncertainty of the experiment as mentioned in Nobel (2013) is 20%, see the red uncertainty bar in Fig. 11. Several conclusions can be drawn from the numerical results: • Both the CFD model and the experimental work show a decreasing trend of dimensionless soil cavity depth as the jet traverse velocity increases. • The trend lines are almost parallel to each other, meaning the slope of the trend lines are close to each other.
The computed width of the cavity w is 0.05 m with an uncertainty of one grid cell size of 0.005 m for all simulations. Nobel (2013) observed that for penetrating jets and 12 < p j /SU < 200 that the cavity width is about 1-1.5 times the nozzle diameter. The dimensionless cavity width in the present simulations is slightly larger ≈ 1.7. Given the experimental uncertainties in e.g. the shear strength and the determination of the cavity dimensions after jetting, the prediction of the cavity depth is reasonably accurate.

Soil surface dislodgement rate
In this section, the soil surface dislodgement rate calculated from CFD simulation will be compared with experimental results of Nobel (2013). CFD Sim 1, which is a penetrating jet case, is used here for validation. Fig. 12 shows the amount of soil dislodged during the experiment. It can be seen that around 9 cells of 5 mm × 5 mm are removed after 0.004 s of jetting. Because of the limitation of the experimental data at hand, the instantaneous surface soil dislodgement rate of the experiment cannot be calculated. With the experimental data at hand, the averaged surface soil dislodgement rate is calculated as 0.056 m 2 /s.
The surface dislodgement rate of the cohesive soil as a function of time calculated by CFD simulation is shown in Fig. 13. The green dash line represents the soil surface dislodgement rate of the experiment. It is found that the soil surface dislodgement rate of the experiment is in between the maximum and minimum value calculated by CFD model. It can also be seen that the dislodgement rate calculated by CFD model has a cyclic behavior with a period of 0.02 s. This relates to the layer-wise failure of the cavity wall as observed in Fig. 12. The jet first needs to step a certain distance on the top of the cavity wall to reach a critical momentum transfer level to initiate shear plane formation inside the wall. After failure the jet has to translate a certain distance before this critical level is reached again. In Section 3.3 this processes will be examined in further detail.
The vertical speed of the soil front u f of the experiment and CFD simulation is also compared, see Fig. 14. The vertical front speeds of the experiment and CFD simulation are both around 12 m/s.

Analysis based on the numerical model
In this section, some analysis based on this CFD model will be carried out on the following aspects: • The pressure exerting on the soil as a function of time • The shear plane change during jetting as a function of time • The soil wall texture (nerves), which shows the soil failure pattern on the soil cavity wall

Pressure on soil surface over time
The cohesive soil fails when the pressure exerting on it exceeds the bearing capacity. The pressure is expected to build up on the soil surface over time until the pressure exerting on the soil is large enough to deform the soil. Fig. 15 shows the computed pressure profile of a fixed point on the soil surface over time. It is seen that when the jet is impinging upon the soil surface, the pressure is building up. When the pressure is large enough to deform the soil, the soil fails in such a way that the flow can enter the cavity with less obstruction (flow more parallel to the cavity surface) and hence the pressure on the cavity surface will decrease. Fig. 16 and Fig. 17 show the maximum pressure on a horizontal line through the soil over time. The horizontal line is on a vertical distance of 0.19 m from the bottom of the computational domain. At this location the jet starts to deflect by more than five degrees indicating the end of the non-deflection and beginning of the deflection zone.
The maximum pressure, as expected, occurs on the intersection between the soil surface and the horizontal line. The soil is represented by the red region in Figs. 16 and 17. It is also observed that the maximum pressure on the soil changes periodically for both the deflecting jet and penetrating jet case with a periodicity of 0.012 s and 0.045 s respectively. Therefore the simulations clearly demonstrate that the failure mechanism is discontinuous: the jet has to step a certain increment forwards before enough force is exerted by the jet flow to the soil before failure occurs. After failure of a soil increment the static pressure decreases again while the jet moves forward to start a new cycle. Fig. 17 also demonstrates that in between the main failure cycles intermediate clay dislodgement events can occur yielding local pressure peaks. In Section 3.3 a simplified failure model is derived based on these observations.
Comparing the maximum pressure in Figs. 15 and 16, it can be seen that the maximum pressure on the soil surface (in Fig. 15) is smaller than that at a lower position (in Fig. 16). This means that the stagnation point of the jet is located somewhat below the original bed level. The stagnation point can be regarded as the location where the jet starts deflecting. It also means that the soil at the surface fails before the arrival of the jet core. This results in a lower pressure at the original bed level. Deeper in the cavity where the jet core hits the clay surface a higher pressure builds up.  (Nobel, 2013). (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)

Transients in shear rate distribution
The shear plane is the demarcation between deforming and nondeforming regions. The exact position of the shear plane cannot be determined from the simulations as the soil is modelled as a regularized Bingham fluid. However, an indication of the location of the shear plane can be obtained by plotting the second invariant I 2 of the shear rate tensor, see Fig. 18. This distinguishes regions with large deformation rates from regions where soil deformation is practically absent. Fig. 18 shows the shear plane patterns of the penetrating jet in Sim 3 over time. The shear plane pattern is marked by the white curve in Fig. 18. It is observed that the shear planes in the penetrating jet regime of Sim 3, change periodically in 0.045 s. This coincides the period of maximum pressure on the soil in Fig. 17. For the deflecting jet in Sim 2E, the shear plane pattern also changes periodically. The period is 0.012 s, which coincides the period of maximum pressure on the soil in Fig. 16.
The failure patterns on the cavity wall for two different experimental conditions (Sim 2E and Sim 3) are shown in Fig. 19. These patterns confirm the layer wise failure of the cavity wall. It can be seen that the patterns on the wall are repeating. This reveals the periodic behavior of the jetting process. The period of the jetting process can be calculated as the pattern distance divided by the jet trail velocity v t . The corresponding periods in Figs. 19a and b are 0.012 s and 0.045 s respectively, which is consistent with the time traces of the maximum pressure in section.3.2.1. The horizontal length scale of texture carved in the cavity wall in Sim 2E and Sim 3 is approximately the same. Sim 3 (Fig. 19b) also demonstrates several pattern repetitions along the jet flow trajectory, which are absent in Sim 2E Fig. 19a. This can be related to the number of sub maxima of the pressure within on cycle (Figs. 16 and 17). This implies that in the penetrating jet regime the cavity depth is reached by a sequence of pressure build up and failure events.

Simplified failure model
To understand the cyclic behaviour of the jet excavation process it is helpful to schematize the problem (Fig. 20). It is assumed that the jet has to travel a certain distance Δx on the soil element (dashed) before the momentum transfer from the jet to the soil element is sufficiently large to overcome the strength of the clay and a shear plane forms at the backside of the soil element. After failure the jet has to step further before a new cycle sets in after Δt time units. The objective is to predict the cavity depth h, the vertical dislodgement velocity v c , Δx and Δt as a function of the jet nozzle exit velocity v j , nozzle radius r, shear strength su = τ y and traverse velocity v t . To solve this problem four equations are required. It is assumed that the soil is incompressible (undrained), which implies that: and h = v c Δt.
The vertical momentum flux of the nozzle F j = πr 2 ρ f v 2 j = 2πr 2 p j is balanced by the normal stress σ on the cavity basis and shear stress τ on the sides of dislodged element. For the sake of simplicity it is assumed that during failure the mobilized shear stress equals the yield stress (undrained shear strength) τ = τ y = su. This simplification is only reasonable if the difference between the undrained shear strength at critical state and peak do not differ considerably and does not depend on the shear rate. For the cohesive soils considered here, this difference is limited and falls within the accuracy range to determine the shear strength with various methods (see Nobel (2013)). The momentum balance becomes: where βF j is the fraction of the momentum flux that acts on the soil element, which is compensated by the shear force, while the fraction (1 − β)F j is the fraction of the momentum flux that acts on the cavity basis.
The total jet power that is exerted on the soil element is 1 2 βF j v j . A fraction α of this power is required to overcome the friction at the shear plane and a fraction (1 − α) is required for the bulk deformation of the soil element. In the deflection zone this bulk deformation consists of large plastic deformations and in the non-deflection zone energy is required to form internal shear planes. The energy balance becomes: Solving the closed set of Equations (9)-(12), yields: It follows that the cavity dimensions and period of the failure process are controlled by the momentum and energy coefficients α and β. These parameters cannot be obtained a priori and require advanced modelling techniques to describe both the momentum distribution of the jet flow inside the cavity and the large soil mechanical deformations of the walls. Table 3 reports the energy and momentum distribution coefficients α and β of the three different simulations computed with the cavity depth h and observed periodicity Δt. The values for α indicates that less than 60% of the jet momentum flux is exerted on the soil element. The values for β indicates that less than 40% of the jet power on the soil element is used for shearing the cavity wall. The values of the coefficients show significant variations between the simulations, demonstrating that the momentum and the energy distribution are a function of the jet traverse velocity, clay strength and jet pressure. The nozzle diameter and stand off distance are not varied in this study. These parameters will probably have an effect on the cavity width w and also on the distribution coefficients α and β. This requires further study.
For simulations 2A-D the width of the cavity determined by the CFD models is approximately the same as for case 2E. The periodicity of the jetting process has not been determined for 2A-D. However, if it assumed that α and β do not depend on the jet traverse velocity it is possible to compute the cavity depth. Fig. 21 shows the computed cavity depths for different jet traverse velocities of simulations 2A-E. It is seen that the simplified failure model provides a reasonable prediction of the cavity depth which indicates that the momentum and energy distribution coefficients do not strongly depend on the jet traverse velocity for this case. This might be a useful observation for finding an optimum between the cavity depth and traverse velocity in terms of production (excavated volume per unit of time). However, a more extended study should be performed to draw firm conclusions in this direction.

Conclusions
The proposed CFD model predicts the cavity depth with reasonable accuracy for a range of jet pressures, undrained shear strengths of the soil and traverse velocities of a single nozzle. The second novelty of this work is, that it clearly confirms that the failure process of a highpressure jet in cohesive soil at sufficiently high traverse velocities is essentially a discontinuous process with a given periodicity. The CFD calculations provide very valuable information regarding the shear and pressure evolution during the cavity formation. This diagnostic information, which has not been obtained by experimental studies so far, is very useful to explain this phenomenon. Based on the analysis of the CFD data a simplified failure model has been developed that provides more basic understanding of the failure process by linking the observed step length, periodicity and cavity depth to the control parameters: soil strength, jet pressure and traverse velocity. This simplified approach is also useful to assist practical calculations. More experimental or CFD data is required to calibrate the coefficients of this simplified failure model. Since the CFD model is generically formulated it can be applied for different nozzle configurations in a wide range of  dredging and offshore engineering equipment such as dragheads, trenchers and drilling tools. This is helpful to improve the design of dredging equipment, optimize operational settings and estimate production.  Fig. 20. Simplification of the failure mechanism of a cohesive soil element by a traversing jet. The shadow part represents the soil element that is deformed and removed by shear during a period Δt, h is the soil cavity depth, w is the soil cavity width and Δx is the length, F jet is the momentum flux from the jet, τ is the shear force on the sides of the element, σ the basal pressure on the cavity wall.
The jet traverses from left to right with velocity v t .

Table 3
Coefficients α and β based on h and Δt of the CFD simulations.