Wave-Induced Seabed Response around a Dumbbell Cofferdam in Non-Homogeneous Anisotropic Seabed

Abstract

Cofferdams are frequently used to assist in the construction of offshore structures that are built on a natural non-homogeneous anisotropic seabed. In this study, a three-dimensional (3D) integrated numerical model consisting of a wave submodel and seabed submodel was adopted to investigate the wave–structure–seabed interaction. Reynolds-Averaged Navier–Stokes (RANS) equations were employed to simulate the wave-induced fluid motion and Biot’s poroelastic theory was adopted to control the wave-induced seabed response. The present model was validated with available laboratory experimental data and previous analytical results. The hydrodynamic process and seabed response around the dumbbell cofferdam are discussed in detail, with particular attention paid to the influence of the depth functions of the permeability $K_i$ and shear modulus $G_j$. Numerical results indicate that to avoid the misestimation of the liquefaction depth, a steady-state analysis should be carried out prior to the transient seabed response analysis to first determine the equilibrium state caused by seabed consolidation. The depth function $G_j$ markedly affects the vertical distribution of the pore pressure and the seabed liquefaction around the dumbbell cofferdam. The depth function $K_i$ has a mild effect on the vertical distribution of the pore pressure within a coarse sand seabed, with the influence concentrated in the range defined by 0.1 times the seabed thickness above and below the embedded depth. The depth function $K_i$ has little effect on seabed liquefaction. In addition, the traditional assumption that treats the seabed parameters as constants may result in the overestimation of the seabed liquefaction depth and the liquefaction area around the cofferdam will be miscalculated if consolidation is not considered. Moreover, parametric studies reveal that the shear modulus at the seabed surface $G_{z0}$ has a significant influence on the vertical distribution of the pore pressure. However, the effect of the permeability at the seabed surface $K_{z0}$ on the vertical distribution of the pore pressure is mainly concentrated on the seabed above the embedded depth in front and to the side of the cofferdam. Furthermore, the amplitude of pore pressure decreases as Poisson’s ratio ${\mu }_s$ increases.

1. Introduction

Nowadays, with the continuous development of ocean engineering, large numbers of offshore structures are being rapidly developed in coastal regions. With the extremely harsh and complex marine environment, long wave periods and high wave energy—among other conditions—major technical challenges arise during the construction of marine structures. In an effort to mitigate such challenges, cofferdams are frequently used to facilitate the construction of offshore structures, such as large-scale sea-spanning bridges [1,2,3,4], because their feasibility and security are superior compared with other available alternatives. In engineering practice, cofferdams are generally manufactured and assembled on the shore and then floated to designated locations, where they are sunk to designated elevations. After sinking, the cofferdam bottom is sealed with concrete with a specified thickness. However, as a large-scale structure, a cofferdam inevitably has a great disturbance effect on the surrounding flow field; in consequence, the seabed response around the cofferdam is significantly affected. Field observations of some project cases have shown that wave loads do not usually cause direct damage to the strength or stiffness of an offshore structure but they often lead to the seabed liquefaction of the region that surrounds the offshore structure and subsequently result in the partial or total loss of bearing capacity [5,6,7,8]. Pile foundations, pipelines and breakwaters have received great attention in the study of wave–seabed–structure interactions [9,10,11], whereas cofferdams have been neglected in such investigations. Therefore, the study of the wave-induced seabed response around cofferdams is of great significance to the safety control during civil engineering construction.

Compared with the traditional uniform seabed, variations in two important soil parameters, permeability and shear modulus, with different burial depths have increasingly become the subject of researchers’ attention. As a consequence of soil consolidation due to both the overburden soil pressure and the water pressure above the water–soil interface, the permeability and shear modulus of soil vary strongly with depth [12,13,14,15,16,17]. For example, the permeability of sediments in the Gulf of Mexico was reported to vary with burial depth [18]. Moreover, the rigidity of soil was found to significantly increase with depth because of the increasing effective overburden pressure [19,20]. Therefore, in order to accurately analyze the seabed response around a cofferdam, the variations in permeability and shear modulus with the buried depth cannot be ignored and the consideration of them will help to correctly estimate the stability of the seabed around the cofferdam.

In recent years, the wave-induced seabed response in the vicinity of offshore structures has been the focus of numerous studies [21,22,23,24,25]. In particular, the wave-induced seabed response around pile foundations [26,27,28], breakwaters [29,30,31] and pipelines [32,33,34,35,36] has attracted a lot of attention. It should be noted that, in the aforementioned studies, the ocean seabed was regarded as a homogeneous porous medium to describe the mechanisms of the wave-induced seabed response around pile foundations, breakwaters and pipelines. Moreover, the soil permeability and shear modulus were assumed to be constant. Thus, inadequate attention has been paid to the depth-related variation in seabed parameters in an inhomogeneous natural seabed. While considering the nonlinearity of the distribution of seabed parameters, Zhou et al. [37] studied the wave-induced seabed response around a pipeline buried in a multilayered seabed and the numerical results showed that pore pressure generation and liquefaction were significantly affected by the permeability of the multilayered seabed. However, the investigations that treated the vertical heterogeneous seabed as a multilayered seabed medium are valid only for uniform materials in each sublayer [38,39,40,41]. In addition, the boundary conditions used at the interface of two sublayers lead to discontinuity in some seabed responses, such as horizontal effective normal stresses, which is inconsistent with the fact that the seabed response varies continuously with depth in a natural seabed. Moreover, treating the seabed as a cross-anisotropic medium, Zhao et al. [42] studied the influences of the cross-anisotropic soil behavior on the wave-induced seabed response around a pipeline in a seabed. In brief, regardless of whether the seabed has been considered uniform, an inhomogeneous multilayered medium or a cross-anisotropic medium, seabed parameters such as permeability and shear modulus have been treated as uniform in previous studies.

To overcome the above deficiencies, Jeng and Lin [14] and Jeng and Seymour [15,16] regarded seabed permeability and/or shear modulus as functions of burial depth and analyzed the variations in wave-induced pore water pressure and effective vertical stress. Their numerical results showed that the soil permeability significantly affected the wave-induced seabed response, especially in a coarser seabed and shear modulus had similarly significant effects on the seabed response in finer seabeds. However, the existence of the structure was taken into account despite the well-known fact that wave propagation and the distribution of the wave-induced seabed response are significantly affected by the presence of a large-scale structure [27]. Therefore, when studying the wave-induced seabed response around a cofferdam, it is necessary to consider the variation in soil parameters in response to the burial depth.

On the basis of the deficiencies reflected by the above-described research situation, in this study, a 3D integrated numerical model was established to investigate the hydrodynamic process of a nonlinear wave–dumbbell-cofferdam interaction and the associated wave-induced seabed response around the dumbbell cofferdam. The wave submodel with the dumbbell cofferdam is governed by Reynolds-Averaged Navier–Stokes (RANS) equations and the seabed submodel is governed by Biot’s poroelastic theory. A detailed description of the integrated numerical model is illustrated in Section 2 and the model’s validation, which was performed using available laboratory experimental data and previous analytical results, is described in Section 3. Then, Section 4 reports the simulations that were carried out using the validated model, namely, the hydrodynamic process of the wave–dumbbell-cofferdam interaction and the wave-induced dynamic seabed response. The particular focus of the parametric study is on the vertical distributions of pore pressure in front of, behind, and to the side of the dumbbell cofferdam.

2. Numerical Model

2.1. Model Description

Figure 1 depicts the wave–seabed–dumbbell-cofferdam interaction. As shown in Figure 1, x, y and z are the Cartesian coordinates. A 3D wave normally propagates to the dumbbell cofferdam from the wave inlet side to the wave outlet side along the positive direction of the x-axis. The integrated model consists of two parts: the wave submodel and seabed submodel. The wave submodel represents the part of the positive region of the z-axis and the seabed submodel represents the part of the negative region of the z-axis. d${}_e$ is the embedded depth of the cofferdam, h is the seabed depth, d is the water depth, L is the wave length and H is the wave height. L${}_x$, L${}_y$ and L${}_z$ are the length, width and height of the computational domain in the x-, y- and z-directions, respectively, and their values should be based on pretest results to minimize the influences of boundary conditions on the flow field around the dumbbell cofferdam.

The dumbbell cofferdam is taken as an impermeable rigid object that is embedded in a porous seabed. A diagram of the dumbbell cofferdam is shown in the local enlargement of Figure 1a, in which locations A–G are seven typical locations at a distance of 1 m from the surface of the dumbbell cofferdam. A, C and E are located in front of the dumbbell cofferdam and B, D and F are behind it, relative to the wave propagation direction; G is on the lateral side of the dumbbell cofferdam at a 90${}^{\circ }$ angle to the wave propagation. The wall thickness is 2 m and concrete is used to seal the bottom of the cofferdam with a thickness of 6 m.

In the wave submodel, a wave-absorbing layer with a certain length is set at the wave outlet to reduce the influence of wave reflection on the flow field. In the seabed submodel, the variable permeability and shear modulus are taken as a function of burial depth because of consolidation under overburden pressure. Three different types are considered, that is, constant, linear and exponential functions, as shown in Figure 2. In a follow-up study, different combinations of variable permeability and shear modulus are considered to assess the influence of variable permeability and shear modulus on the wave-induced seabed response around the dumbbell cofferdam. Type${}_{ij}$ denotes the combination of Type K${}_i$ and G${}_j$. It is important to note that Type${}_{11}$ is the conventional assumption for a homogeneous seabed. Moreover, in this study, the value of permeability at the rigid bottom (z = −h) is assumed to be one-tenth of that at the seabed surface K${}_{z0}$ and the value of the shear modulus at the rigid bottom is assumed to be ten times larger than that at the seabed surface G${}_{z0}$ [14]. The parameters used in the following numerical experiments are listed in

Table 1. Input data for the standard numerical simulation.

Parameters	Value	Unit
Wave parameters
Wave height (H)	5.0	m
Water depth (d)	25.0	m
Wave period (T)	10.0	s
Seabed parameters
Density of soil (${\rho }_s$)	2650	kg/m${}^3$
Permeability at seabed surface ($K_{z0}$)	1 × 10${}^{-2}$, 1 × 10${}^{-4}$	m/s
Shear modulus at seabed surface ($G_{z0}$)	5.0 × 10${}^6$	N/m${}^2$
Poisson’s Ratio (${\mu }_s$)	0.33	-
soil porosity ($n_s$)	0.4	-
Degree of saturation ($S_r$)	0.98	-
Seabed thickness (h)	50.0	m
Dumbbell cofferdam parameters
Cofferdam length ($l_c$)	60.0	m
Embedded depth ($d_e$)	15.0	m
Cofferdam density (${\rho }_c$)	2400	kg/m${}^3$
Shear modulus ($G_c$)	1.25 × 10${}^{10}$	N/m${}^2$
Poisson’s Ratio (${\mu }_c$)	0.2	-

.

2.2. Governing Equations

2.2.1. Wave Submodel

Wave propagation is modeled by the RANS equation and the $\kappa$-$ϵ$ turbulence model. As a three-dimensional model, the incompressible fluid motion can be described by the mass equation and the RANS equation with a $\kappa$-$ϵ$ closure [43,44], which can be expressed as

$$\frac{\partial 〈u_i〉}{\partial x_i}=0,$$

(1)

$$\frac{\partial \rho 〈u_i〉}{\partial t}+\frac{\partial \rho 〈u_i{u}_j〉}{\partial x_j}=-\frac{\partial }{\partial x_i}\left[〈p〉+\frac{2}{3}\rho \kappa \right]+\frac{\partial }{\partial x_j}\left[{\mu }_{eff}\left(\frac{\partial 〈u_i〉}{\partial x_j}+\frac{\partial 〈u_j〉}{\partial x_i}\right)\right]+\rho g_i,$$

(2)

in which $x_i$ is the Cartesian coordinate, $〈u_i〉$ is the ensemble mean velocity component, t is the time, $\rho$ is the fluid density, $〈p〉$ is the fluid pressure and $g_i$ is the acceleration of gravity. ${\mu }_{eff}=\mu +{\mu }_t$ is the total effective viscosity, $\mu$ is the dynamic viscosity, ${\mu }_t$ is the turbulent viscosity and $\kappa$ is the turbulence kinetic energy.

The volume of fluid (VOF) method [45] is applied in this model to track the free water surface. The VOF method aims to define a function F to represent the fractional volume of water fluid:

$$\frac{\partial F}{\partial t}+\frac{\partial \left(u_{i}F\right)}{\partial x_i}=0,$$

(3)

In Equation (3), $F=1$ indicates that the cell is full of water, while $F=0$ corresponds to a cell fully occupied by air. A cell with a value of $0<F<1$ contains a free water surface.

2.2.2. Seabed Submodel

In the seabed submodel, the porous seabed is considered to be elastic, unsaturated and hydraulically permeable. Moreover, the soil skeleton and the pore fluid are assumed to be compressible and to obey Hooke’s law. Biot’s poroelastic theory [46] is used to govern the soil response, which can be expressed as

$$\frac{K_x}{K_z}\frac{{\partial }^2{p}_s}{\partial x^2}+\frac{K_y}{K_z}\frac{{\partial }^2{p}_s}{\partial y^2}+\frac{{\partial }^2{p}_s}{\partial z^2}+\frac{1}{K_z}\frac{d{K}_z}{dz}\frac{\partial p_s}{\partial z}-\frac{{\gamma }_w{n}_s{\beta }_s}{K_z}\frac{\partial p_s}{\partial t}=\frac{{\gamma }_w}{K_z}\frac{\partial {ϵ}_s}{\partial t},$$

(4)

where $p_s$ is the wave-induced pore pressure, ${\gamma }_w$ is the unit weight of water, $n_s$ is the soil porosity and $K_x$, $K_y$ and $K_z$ are the soil permeability in the x-, y- and z-directions, respectively. Note that the wave-induced pore pressure is insensitive to the hydraulic anisotropic parameters $K_x/K_z$ and $K_y/K_z$ over the range of 1–100 [14]; hence, only isotropic conditions are considered. ${\beta }_s$ is the compressibility of the pore fluid and ${ϵ}_s$ is the volume strain, which can be expressed as

$${\beta }_s=\frac{1}{K_w}+\frac{1-S_r}{P_{w0}},\mathrm{and}{ϵ}_s=\frac{\partial u_s}{\partial x}+\frac{\partial v_s}{\partial y}+\frac{\partial w_s}{\partial z},$$

(5)

where $K_w$ is the bulk modulus of pore water, $S_r$ is the degree of seabed saturation and $P_{w0}$ is the absolute static pressure. $u_s$, $v_s$ and $w_s$ are the soil displacements in the x-, y- and z-directions, respectively.

The shear modulus of the soil $G\left(z\right)$ is related to Young’s modulus $E\left(z\right)$ by Poisson’s ratio ${\mu }_s$ in the form of $E\left(z\right)$/2(1 + ${\mu }_s$). Since $G\left(z\right)$ is assumed to be a function of burial depth z here, according to linear elasticity theory, the equations of force equilibrium can be expressed as

$$G\left(z\right){\nabla }^2{u}_s+\frac{G\left(z\right)}{1-2{\mu }_s}\frac{\partial {ϵ}_s}{\partial x}+\frac{dG\left(z\right)}{dz}\left[\frac{\partial u_s}{\partial z}+\frac{\partial w_s}{\partial x}\right]=-\frac{\partial p_s}{\partial x},$$

(6)

$$G\left(z\right){\nabla }^2{v}_s+\frac{G\left(z\right)}{1-2{\mu }_s}\frac{\partial {ϵ}_s}{\partial y}+\frac{dG\left(z\right)}{dz}\left[\frac{\partial w_s}{\partial y}+\frac{\partial v_s}{\partial z}\right]=-\frac{\partial p_s}{\partial y},$$

(7)

$$G\left(z\right){\nabla }^2{w}_s+\frac{G\left(z\right)}{1-2{\mu }_s}\frac{\partial {ϵ}_s}{\partial z}+\frac{dG\left(z\right)}{dz}\left[\frac{{\mu }_s}{1-2{\mu }_s}\left(\frac{\partial u_s}{\partial x}+\frac{\partial v_s}{\partial y}\right)+\frac{2(1-{\mu }_s)}{1-2{\mu }_s}\frac{\partial w_s}{\partial z}\right]=-\frac{\partial p_s}{\partial z},$$

(8)

in the x-, y-, z-directions, respectively.

Moreover, the modified Biot’s theory is used to describe the motion of the dumbbell cofferdam (neglecting the pore pressure) and can be expressed as

$$G_c{\nabla }^2{u}_c+\frac{G_c}{1-2{\mu }_c}\frac{\partial {ϵ}_c}{\partial x}=0,$$

(9)

$$G_c{\nabla }^2{v}_c+\frac{G_c}{1-2{\mu }_c}\frac{\partial {ϵ}_c}{\partial y}=0,$$

(10)

$$G_c{\nabla }^2{w}_c+\frac{G_c}{1-2{\mu }_c}\frac{\partial {ϵ}_c}{\partial z}=0,$$

(11)

where $u_c$, $v_c$ and $w_c$ are the displacements of the cofferdam in the x-, y-, z-directions, respectively. ${ϵ}_c$ = $\partial u_c$/$\partial x$ + $\partial v_c$/$\partial y$ + $\partial w_c$/$\partial z$ is the volume strain of the cofferdam. $G_c$ is the shear modulus and ${\mu }_c$ is the Poisson’s ratio of the cofferdam.

2.3. Boundary Conditions

The oscillatory pore pressure and the displacements within the soil can be obtained by solving the governing equations introduced above with approximate boundary conditions. To calculate the water wave pressure acting on the seabed surface, relevant boundary conditions need to be specified in the wave submodel. Similarly, the boundary conditions should also be specified in the seabed submodel when determining the wave-induced seabed response around the dumbbell cofferdam.

2.3.1. Wave Submodel

In the wave submodel, appropriate boundary conditions are required to solve the RANS equations to obtain the water wave pressure $p_w$, which is needed in the seabed submodel:

(1)

The waves are generated at the wave inlet boundary at the position of $x=0$. The outflow boundary is applied at the wave outlet area at $x=L_x$ and a sponge layer is adopted to absorb waves to minimize the wave reflection, which can induce the rising or falling of the free water surface during wave propagation simulation.

(2)

For the two lateral boundaries ($y=0$ and $y=L_y$), a zero-gradient condition (symmetric boundary) is applied, which means that the wave can slide freely along the boundary without penetration.

(3)

At the bottom of the numerical wave tank ($z=0$), the wall boundary is adopted to ensure that the velocity of the fluid in the normal direction is zero. At the vertical maximum position of the computational domain ($z=L_z$), the static pressure condition is set by using fluid elevation.

(4)

At the free water surface, the pressure is set to be equal to the atmospheric pressure, 101 kPa. At the wave–structure boundary, a nonslip condition for velocity is imposed.

2.3.2. Seabed Submodel

(1)

At the seabed surface, the vertical effective normal stress ${\sigma }_z^{\prime }$ and the shear stress vanish and the dynamic pore pressure is equal to the dynamic wave pressure [47]:

$${\sigma }_z^{\prime }={\tau }_{xz}={\tau }_{yz}=0\mathrm{and}{p}_s=p_w-{\gamma }_{w}d\mathrm{at}z=0,$$

(12)

where wave pressure $p_w$ can be obtained from the wave submodel.

(2)

At the bottom of the seabed, an impermeable rigid boundary condition is applied in which the seabed displacements are zero and no vertical flow occurs. This can be expressed as

$$u_s=v_s=w_s=0\mathrm{and}\frac{\partial p_s}{\partial z}=0\mathrm{at}z=-h.$$

(13)

(3)

At the four lateral boundaries of the seabed, the materials are assumed to be impermeable and rigid. Therefore, the horizontal displacements of the four lateral boundaries are set at zero and the normal gradient of the pore-water pressure is zero, which can be expressed as

$$u_s=v_s=0,\frac{\partial p_s}{\partial n}=0.$$

(14)

(4)

At the seabed–structure boundary, the normal gradient of the pore pressure is zero at the impermeable and rigid structure surface. Moreover, there is no relative displacement of soil with respect to the structure, as expressed by

$$\frac{\partial p_s}{\partial n}=0,u_{soil}=u_{structure},v_{soil}=v_{structure},\mathrm{and}{w}_{soil}=w_{structure}.$$

(15)

2.4. Integration of the Wave and Seabed Submodels

In this study, an integrated model that comprises a wave submodel and seabed submodel was constructed. The wave submodel simulates the wave motion and obtains the dynamic wave pressure acting on the porous seabed and the dumbbell cofferdam at different moments. The CFD code FLOW-3D was applied to establish the numerical wave tank and the boundary velocity inlet method was applied to generate waves. The seabed submodel was established by the finite element software COMSOL Multiphysics. The Biot quasi-static governing equations of the seabed submodel can be input through the interface of the PDE module. Once the wave pressure on the seabed surface is calculated from the wave submodel, it is taken as the boundary condition of the seabed submodel to determine the soil response. On the basis of the updated boundary conditions and initial value conditions, the seabed submodel can obtain the pore pressure, displacement and effective stresses.

In the wave submodel, a structured mesh is defined in a Cartesian Coordinate System, as shown in Figure 3a. The dumbbell cofferdam is embedded in the mesh by partially blocking cell volumes and face areas to make the mesh and geometry independent. Moreover, local refinement of the grid around the cofferdam captures the complex flow domains around the dumbbell cofferdam. The fine mesh block contains about 19.98 million cells in total. The length of the wave submodel is 7 times the wave length L and comprises the 6 L length of the wave flume and 1 L length of the wave-absorbing layer (Figure 1), which is used to eliminate the effects of the boundary conditions and wave reflection on the flow field. Moreover, in the seabed submodel, shown in Figure 3b, unstructured Lagrange elements with a maximum global size of 12.5 m in the region away from structure and a minimum size of 0.6 m in the region near the structure are adopted to discretize the whole seabed submodel, including the seabed soil and dumbbell cofferdam. The refinement of the mesh is carried out around the dumbbell cofferdam and the seabed submodel contains 167,779 cells in total. The length of seabed model is about 6 L, which has been assessed as having sufficient accuracy for analyzing the seabed response around a dumbbell cofferdam [48].

3. Model Validation

In this section, the validation of the numerical model is presented. The model was validated by comparing it with previous laboratory experimental data and numerical results. The verification content includes the following three aspects:

(1)

Free surface elevation and wave loading in the wave–structure interaction to verify the ability of the present numerical model to simulate turbulence around the structure.

(2)

Wave-induced seabed response of a non-homogeneous anisotropic seabed to verify the necessity of considering variable permeability and shear modulus when studying the seabed response.

(3)

Seabed oscillatory pore pressure of a non-homogeneous anisotropic seabed to verify the ability of the present numerical model to determine the oscillatory response of a non-homogeneous anisotropic seabed.

3.1. Free Surface Elevation and Wave Loading

Mo et al. [49] studied the interaction between periodic waves and a vertical slender pile in a large wave flume. The test setup is shown in Figure 4. To measure water surface elevation, four wave-gauges were located around the circular cylinder and one was arranged near to the wave flume wall in the direction of wave propagation. The comparisons are presented here for a water depth d = 4.76 m, a wave period T = 4 s, a wave height H = 1.2 m and a cylinder diameter D = 0.7 m.

Figure 5 shows the surface elevations measured at four wave-gauges. The surface elevations are non-dimensionalized by the maximum elevation of $W{G}_{ref}$ and denoted by ${\eta }_{max,W{G}_{ref}}$. Figure 5 exhibits a reasonably good agreement between the measured data and the present numerical results. There is even a good prediction of some characteristic ripples, which were measured in the front of the cylinder (named $W{G}_1$) while the wave trough passes the pile.

Figure 6 shows the time-history curves of the wave loading acting on the structure and the wave loading $F_x$ is non-dimensionalized by $\rho gD{H}^2$. The results show that the numerical results of wave loading are in good agreement with each other in terms of amplitude and period.

In general, the present numerical model can well simulate the turbulent flow around the structure and the wave loading impact on the structure, which in turn can provide reliable data to determine the subsequent wave–structure interaction.

3.2. Wave-Induced Seabed Response in Non-Homogeneous Anisotropic Seabed

Jeng and Lin [14] concentrated on the short-crested-induced soil response in a porous seabed and took the variable permeability and shear modulus as a function of burial depth, as shown in Figure 2. A total of five different combinations of variable permeability and shear modulus, that is, $Typ{e}_{i1}$ and $Typ{e}_{1j}$ (i or j = 1, 2, 3), are considered to verify the necessity and importance of considering a non-homogeneous seabed. It is important to note that $Typ{e}_{11}$ is the conventional assumption for a homogeneous seabed.

Figure 7 is the vertical distribution of the maximum pore pressure in a non-homogeneous anisotropic seabed. The pore pressure is non-dimensionalized by the maximum pore pressure $P_0$ impact on the seabed surface on the basis of linear theory. The results show that the present numerical results are in good agreement with Jeng’s finite element solutions for the numerical experiments, except for some unavoidable deviation caused by the different algorithms used. Moreover, the comparison results show that the present numerical model is accurate and reliable for studying the influence of variable permeability and shear modulus on the seabed response. On the other hand, the soil permeability and shear modulus can significantly affect the wave-induced seabed response.

3.3. Oscillatory Pore Pressure of Non-Homogeneous Anisotropic Seabed in the Experiment

The oscillatory pore pressure of a non-homogeneous anisotropic seabed is compared with the experimental data presented by Liu [50] to verify the importance of considering variable permeability and shear modulus. The seabed consists of 20 sublayers and the soil parameters (permeability coefficient and shear modulus) measured in the experiment are shown in Figure 8. The comparisons are presented here for a water depth d = 5.2 m, a wave period T = 9 s and a wave height H = 5.2 m.

Figure 9 shows the comparisons of the measured and predicted vertical distribution of the maximum pore pressure. As seen in the figure, the numerical results are in good agreement with the experimental results except for minor deviations, which can be attributed to the seabed settlement in the experimental process. During the experiment, the seabed parameters of the non-homogeneous anisotropic seabed, such as permeability and shear modulus, are ongoing changes with the experimental process, which are regarded as constants in the numerical simulation without change.

In general, the present numerical model is capable of further simulating the wave–layered-seabed interaction.

4. Results and Discussion

In this section, the results of the simulations using the validated model are described to illustrate the hydrodynamic process and seabed response resulting from wave action, with a primary focus on the influence of variable permeability and shear modulus. The parameters adopted in the numerical experiments are listed in Table 1 unless specified otherwise. Seven typical locations were selected around the dumbbell cofferdam, as shown in Figure 1a, to more accurately study the strong nonlinearity around the dumbbell cofferdam during interactions between the wave and dumbbell cofferdam.

4.1. Hydrodynamic Process Involved in the Wave and Dumbbell Cofferdam Interactions

To study the process of the wave–dumbbell-cofferdam interaction, a more accurate simulation of the flow field around the dumbbell cofferdam and wave force impact on the dumbbell cofferdam is needed.

Figure 10 shows snapshots of wave propagation and interactions with the dumbbell cofferdam. Four typical times are considered, namely, 78.5, 80.0, 81.0 and 83.5 s. The wave crest and wave trough propagate to the front of the dumbbell cofferdam at t = 78.5 and 83.5 s, respectively. The graphical results show that the wave surfaces at the position of the wave-absorbing layer near the outlet boundary of the numerical flume are approximate zero, which indicates that the wave-absorbing layer has a good wave dissipation effect and can greatly reduce the influence of wave reflection on the calculation results during a long-time numerical analysis. Moreover, as reflected clearly on the contour lines, the wave surface in front of the dumbbell cofferdam expands or shrinks locally and remarkably in the grab handle of the dumbbell cofferdam. Around the dumbbell cofferdam, there emerge pronounced diffraction phenomena caused by the masking effect of the large dumbbell cofferdam on the incident waves.

Figure 11 shows the time histories of the wave profiles at seven typical locations (Locations A–G) around the dumbbell cofferdam. The results show that the existence of the dumbbell cofferdam has a notable influence on the flow field. The masking effect of the large-sized dumbbell cofferdam results in a phase difference in the wave surface at the typical locations and the amplitude of the wave profiles in front of the structure is larger than that on the side and in the back of the dumbbell cofferdam; this observation is consistent with the results shown in Figure 10 and further reflects the superposition process of the incident wave, reflection wave and diffraction wave during the wave–dumbbell-cofferdam interaction. Moreover, the flow field is stable in the time range of 78.5–83.5 s, which was selected as the research time.

Besides the concern about the wave profiles around the dumbbell cofferdam, accurately modeling the wave loading impact on the dumbbell cofferdam is crucial. Figure 12 is the dynamic pressure impact on the dumbbell cofferdam in the wave propagation direction. As shown in the figure, the action area of the dumbbell cofferdam varies with the wave propagation; this variation is one of the main reasons for the change in wave force impact on the dumbbell cofferdam. Moreover, the amplitude of the wave force exerted on the dumbbell cofferdam is 76,117 kN.

The aforementioned simulated results indicate that the wave field around the dumbbell cofferdam is highly three-dimensional and obviously nonlinear and the wave loading impact on the dumbbell cofferdam is significantly related to the wave fluctuation around the dumbbell cofferdam. Moreover, the existence of the dumbbell cofferdam can modify the propagation pattern of its nearby wave flow, which will subsequently affect the dynamic soil response in the seabed.

4.2. Consolidation of the Seabed

The construction process of the dumbbell cofferdam inevitably disturbs the surrounding seabed, which has completed consolidation under gravity and reaches a new equilibrium state. Therefore, it is necessary to first carry out steady-state analysis in the numerical model to take into account the initial consolidation state of the surrounding seabed under the gravitational action of the cofferdam [22,51]. Then, steady-state analysis is used as the initial condition for further analysis of the wave-induced non-homogeneous anisotropic seabed response around the cofferdam.

Figure 13 shows the displacement distributions of the seabed around the cofferdam after the consolidation of the non-homogeneous anisotropic coarse sand seabed for $Typ{e}_{33}$. The results show that seabed displacements are largely concentrated around the dumbbell cofferdam and that the vertical displacement is dominant, as shown in Figure 13a. Because of the large vertical displacement, the lateral displacements near the seabed surface move in different directions relative to the cofferdam and in directions opposite of those of the lateral displacements near the bottom edge of the dumbbell cofferdam, as shown in Figure 13b. Moreover, because of the gravity of the cofferdam, the soil beneath and around the cofferdam is compressed and the initial effective stress beneath and around the cofferdam increases, as shown in Figure 13c. Moreover, steady-state analysis of consolidation can prevent the misestimation of liquefaction development of the soil beneath the structures, considering the local variation in the effective stress around the dumbbell cofferdam. These variations are discussed in detail in the following section.

Figure 14 describes the maximum vertical displacements of the different seabeds (denoted as $K_i$ and $G_j$) after consolidation due to the cofferdam’s self-weight, normalized with the maximum vertical displacement of $Typ{e}_{11}$. The influence of the depth function of permeability $K_i$ on the vertical displacement caused by the consolidation process can be neglected. However, the vertical displacement due to consolidation is more sensitive to the change in shear modulus and the larger the shear modulus, the smaller the vertical displacement.

4.3. Seabed Pore Pressure within the Non-Homogeneous Anisotropic Seabed around the Dumbbell Cofferdam

The hydrodynamic process involved in the wave and dumbbell cofferdam interactions imposes periodically changing dynamic pressures on the seabed surface and this causes the subsequent variation in the seabed response. Computational results are presented here for the oscillatory pore pressure. The typical time at which the wave trough impacts the cofferdam front was selected as t = 83.5 s for further study, given that seabed stability is weak at this moment. In this section, the discussion of the numerical study is concentrated on the seabed pore pressure at the characteristic locations (shown in Figure 1a, Locations A–G), as shown in Figure 15. It can be seen from Figure 15 that the spatial variations in pore pressure are concentrated in the seabed layer above the embedded depth of the cofferdam and the variations are small below the embedded depth for the coarse sand seabed. Moreover, the pore pressure variations at the characteristic locations in front of the dumbbell cofferdam (Locations A, C and E) are more pronounced those in behind it (Locations B, D and F). This is a result of the blocking effect of the dumbbell cofferdam on the wave trough propagating forward; the blocking effect decreases the wave surface near the handle of the dumbbell cofferdam and increases the suction effect. Because of the similarity in pore pressure distributions at locations A and C, the characteristic locations A, B, E, F and G were selected to study the vertical distributions of pore pressure in front of, behind and at the lateral side of the dumbbell cofferdam and to then study the influence of a non-homogeneous anisotropic seabed on the vertical distributions of pore pressure.

Figure 16 and Figure 17 depict the vertical distributions of the pore pressure at characteristic locations at t = 83.5 s. The pore pressure $p_s$ is non-dimensionalized by the amplitude of the wave pressure $P_0$ calculated by linear wave theory. It is seen that the depth function of shear modulus $G_j$ obviously affects the vertical distributions of the pore pressure: with the increase in shear modulus, the amplitude of the pore pressure decreases in the same location. However, the depth function of permeability $K_i$ has a mild effect on the coarse sand seabed and the influence concentrates near the embedded depth; however, it has little effect on a fine sand seabed. Moreover, two points of sudden change in the vertical distributions of pore pressure emerge: one point near the seabed surface ($z/h$ = −0.05) and another point near the cofferdam embedded depth ($z/h$ = −0.3). For the coarse sand seabed, shown in Figure 16, the first point ($z/h$ = −0.05) at the characteristic locations behind the dumbbell cofferdam (Locations B and F) is more obvious than that at the characteristic locations in front of the dumbbell cofferdam (Locations A and E). At the second point ($z/h$ = −0.3), the effect of the depth function $K_i$ on the vertical distribution of pore pressure in the range of 0.1h above and below the embedded depth is the opposite. In the range between the embedded depth and 0.1h above it, with the decrease in permeability, the amplitude of the pore pressure decreases (see the local enlargements in Figure 16 for details). For the fine sand seabed, shown in Figure 17, there are some sawtooth phenomena in the vertical distributions of the pore pressure due to the displacement coordination at the interface between the structure and seabed. Moreover, two points of sudden change in the vertical distributions of pore pressure are obvious.

4.4. Seabed Liquefaction within the Non-Homogeneous Anisotropic Seabed around the Dumbbell Cofferdam

It must be noted that seabed liquefaction in the vicinity of a dumbbell cofferdam is a primary concern in offshore engineering practice. In this section, the wave-induced liquefaction potential around the dumbbell cofferdam in a non-homogeneous anisotropic seabed is further investigated. The theory proposed by Jeng [9] considers the spatial effects of the cofferdam, which can be expressed as

$$-\frac{1}{3}\left({\gamma }_s-{\gamma }_w\right)\left(1+2K_0\right)z+\frac{1}{3}\left(\Delta {\sigma }_{x0}^{\prime }+\Delta {\sigma }_{y0}^{\prime }+\Delta {\sigma }_{z0}^{\prime }\right)\le p_s-p_b,$$

(16)

where the left-hand side represents the average effective geostatic stress considering the seabed consolidation due to the cofferdam weight and the right-hand side ($p_s-p_b$) denotes the excess pore pressure. $p_b$ is the dynamic wave pressure and ${\gamma }_s$ and ${\gamma }_w$ are the unit weights of soil and water. $K_0$ is the coefficient of lateral earth pressure, with $K_0$ = 5 being the commonly used value. $\Delta {\sigma }_{x0}^{\prime }$, $\Delta {\sigma }_{y0}^{\prime }$ and $\Delta {\sigma }_{z0}^{\prime }$ are the effective stresses induced by seabed consolidation.

It is noted that the liquefied region occurs under wave troughs rather than wave crests because the wave-induced pore pressure is pushing the soil particle upward near the wave troughs. The moment that the wave trough propagates to the front of dumbbell cofferdam is of interest here. Figure 18 depicts the liquefaction zone with and without considering consolidation (x–z plane, y = 200; y–z plane, x = 392) within the coarse seabed for $Typ{e}_{33}$ at t = 83.5 s. When considering seabed consolidation, it is seen that the maximum liquefaction depth is obtained in front of the dumbbell cofferdam and no liquefaction occurs behind the dumbbell cofferdam, as shown in Figure 18a,b. The liquefied seabed flows as a fluid under the periodic action of waves and because of the difference in liquefaction depth between the front and back of the dumbbell cofferdam on the x–z plane, the dumbbell cofferdam tends to be inclined and to ultimately fail. Without considering seabed consolidation, the maximum liquefaction depth is located on the surface of the cofferdam and the value is larger. On the other hand, in addition to occurring in front of the cofferdam, liquefaction areas also appear both at the bottom of and behind the cofferdam. This does not conform to the actual situation and the deviation is due to the neglect of the redistribution of internal forces around structures as a result of seabed consolidation. Obviously, the liquefaction area around the cofferdam will be miscalculated if seabed consolidation is not considered.

To further study the influence of non-uniformity of the seabed on the spatial distributions of seabed liquefaction around the cofferdam, Figure 19 shows the liquefied zone on the x–z plane (y = 200) at t = 83.5 s. It can be seen that the depth functions of shear modulus $G_j$ significantly influence seabed liquefaction, while the depth functions of permeability $K_i$ have little effect, regardless of whether coarse sand or fine sand was used in the numerical experiments presented here. This is attributed to the weak influence of $K_i$ on the pore pressure within the seabed surface. Moreover, with the increase in the shear modulus, the liquefaction depth around the cofferdam is smaller; for example, the liquefied zone of $Typ{e}_{11}$ is larger than that of $Typ{e}_{13}$. It is also noticeable that a small area that ranges from 0.01 h to 0.05 h near the cofferdam surface is not liquefied. This discontinuity might result from the displacement coordination at the interface between the structure and the seabed, as clearly displayed in the local enlargement in Figure 19b. Therefore, the traditional assumption that treats seabed parameters as constants may overestimate the seabed liquefaction depth.

4.5. Parametric Studies on Pore Pressure Response around the Dumbbell Cofferdam

The aforementioned content focuses on the influence of the depth functions $K_i$ and $G_j$ on the seabed response around the dumbbell cofferdam. In this section, the effect of permeability and shear modulus at the seabed surface, that is, $K_{z0}$ and $G_{z0}$, as well as the Poisson’s ratio ${\mu }_s$, on the vertical distributions of pore pressure at characteristic locations is studied with the same distribution function (taking $Typ{e}_{33}$ as an example). The characteristic locations A, B and G were selected to study the vertical distributions of pore pressure in front of, behind and to the side of the cofferdam, respectively.

Figure 20 shows the vertical distributions of pore pressure with various permeabilities at the seabed surface $K_{z0}$ (10${}^{-1}$, 10${}^{-2}$, 10${}^{-3}$ and 10${}^{-4}$ m/s) for $Typ{e}_{33}$ at t = 83.5 s. Compared with the vertical distributions of pore pressure behind the cofferdam, the effect of $K_{z0}$ on the vertical distribution of pore pressure is mainly concentrated on the seabed above the embedded depth in front of and to the side of the cofferdam. The two above-mentioned points of sudden change at $z/h$ = −0.05 and $z/h$ = −0.3 in the vertical distributions of pore pressure manifest clearly when $K_{z0}$ is relatively small (10${}^{-3}$ m/s for the present numerical experiments). Then, Figure 21 depicts the vertical distributions of pore pressure with various shear moduli at the surface of a coarse seabed $G_{z0}$ (1.0 × 10${}^6$ pa, 5.0 × 10${}^6$ pa, 1.0 × 10${}^7$ pa) for $Typ{e}_{33}$ at t = 83.5 s. It shows that $G_{z0}$ has a significant influence on the vertical distributions of pore pressure. With the increase in $G_{z0}$, the amplitude of the pore pressure is smaller. Further, the two points at $z/h$ = −0.05 and $z/h$ = −0.3 are relatively obvious when the shear modulus is small (1.0 × 10${}^6$ pa for the present numerical experiments). Therefore, the point of sudden change at $z/h$ = −0.3, which reflects the influence of the dumbbell cofferdam on the seabed response, is more significant when $K_{z0}$ and $G_{z0}$ have smaller values. Moreover, The vertical distributions of pore pressure with various Poisson’s ratio ${\mu }_s$(0.33, 0.40, 0.45) are presented in Figure 22. Poisson’s ratio reflects the elastic constant of soil transverse deformation. The results indicate that ${\mu }_s$ has an effect on the vertical distributions of pore pressure. Generally speaking, the amplitude of pore pressure decreases as ${\mu }_s$ increases.

5. Conclusions

In this study, a 3D integrated numerical model was established to investigate the hydrodynamic process of the nonlinear wave–dumbbell-cofferdam interaction and the associated seabed response under wave dynamic loading. On the basis of the numerical results, the following conclusions can be drawn:

(1)

The 3D integrated numerical model is capable of reproducing the 3D hydrodynamic process involved in the wave–structure interaction and wave-induced seabed response. The model was validated using analytical solutions and experimental data and a fairly good agreement was obtained.

(2)

Steady-state analysis should be carried out and then imported as the initial condition for the transient analysis of the seabed response to account for the initial consolidation state due to the dumbbell cofferdam self-weight. On one hand, the initial effective stress beneath and around the cofferdam increases and the steady-state analysis of consolidation can prevent the misestimation of liquefaction development of the soil beneath structures. On the other hand, compared with the depth function of permeability $K_i$, the depth function of shear modulus $G_j$ markedly affects the amplitude of the vertical displacement. This highlights the necessity of considering the vertical distributions of the permeability and shear modulus.

(3)

The seabed response around the dumbbell cofferdam is affected by the combination of the depth functions $K_i$ and $G_j$, that is, $Typ{e}_{ij}$, to different degrees. The spatial variations in pore pressure are concentrated in the seabed layer above the embedded depth of the cofferdam. The depth function $G_j$ obviously affects the vertical distribution of the pore pressure: with the increase in $G_j$, the amplitude of the pore pressure decreases in the same location. However, the depth function $K_i$ has a mild effect on the coarse sand seabed and the influence concentrates in the range defined by 0.1h above and below the embedded depth; however, it has little effect on a fine sand seabed. Moreover, two points of sudden change ($z/h$ = −0.05 and $z/h$ = −0.3) emerge in the vertical distributions of pore pressure and the distinct degree of the two points is related to $Typ{e}_{ij}$.

(4)

Under wave troughs, the liquefaction depth is larger in front of the dumbbell cofferdam than behind it, so the dumbbell cofferdam may tend to be inclined and to ultimately fail as a result of the asymmetry of the liquefaction depth. The depth function $G_j$ has a significant influence on the seabed liquefaction, while the depth function $K_i$ has little effect, regardless of whether coarse sand or fine sand is used in the present numerical experiments. Moreover, the traditional assumption that treats seabed parameters as constants may overestimate the seabed liquefaction depth.

(5)

Parametric studies reveal that the effect of permeability at the seabed surface $K_{z0}$ on the vertical distributions of pore pressure is mainly concentrated on the seabed above the embedded depth in front and to the side of the cofferdam. The shear modulus at the seabed surface $G_{z0}$ has a significant influence on the vertical distribution of pore pressure. With the increase in $G_{z0}$, the amplitude of the pore pressure is smaller. The point of sudden change at $z/h$ = −0.3, which reflects the influence of the dumbbell cofferdam on the seabed response, is more significant when $K_{z0}$ and $G_{z0}$ have smaller values. Moreover, the amplitude of the pore pressure decreases as Poisson’s ratio ${\mu }_s$ increases.