The effects of slab geometries and wave directions on the steep wave-induced soil response and liquefaction around gravity-based offshore foundations

ABSTRACT

Gravity-based offshore foundations generally consist of a bottom slab and one or more cylindrical shafts on top of it. The slab geometry can strongly affect the surrounding flow pattern, dynamic wave pressure distribution and further momentary soil liquefaction. In the present study, the effects of the slab geometry and the incoming wave angle to the foundations on the surrounding soil response are investigated. Gravity-based foundations with bottom slabs of cylindrical shape and hexagonal prismatic shape are considered. For the hexagonal foundation, two different incoming wave angles, i.e. propagating to the hexagon corner and the edge, are simulated. The waves are fully nonlinear steep non-breaking waves. The wave-structure interaction is validated against an existing experiment. Soil consolidation behaviour underneath the foundations is investigated. It is found that both the slab geometry and the incoming wave angle can affect the distributions of pore pressure and momentary liquefaction depth in the surrounding soil.

KEYWORDS:

• Steep wave
• soil response
• momentary liquefaction
• gravity-based
• offshore foundations

1. Introduction

Offshore foundations are installed in an environment that may encounter steep and nonlinear waves. The wave-induced soil liquefaction around the offshore foundations may cause severe soil and further structure failures. The prediction of the wave-induced soil response and liquefaction risk around the offshore foundation is critical in the engineering design. In the recent years, most of the wave-structure-seabed interaction studies have focussed on the coastal and offshore structures such as breakwaters (Jeng 1997; Ulker et al. 2010; Ye 2012b; Celli et al. 2019) and monopile foundations (Li et al. 2011; Chang and Jeng 2014; Zhang et al. 2015; Sui et al. 2016; Lin et al. 2017; Zhao et al. 2017; Zhu et al. 2018), but very few literature studies have considered the wave-structure-seabed interaction of gravity-based foundations with more complex three-dimensional geometries. Recently, Li et al. (2018) studied the seabed consolidation and liquefaction around a hexagonal foundation using an integrated wave-structure-seabed interaction model developed in OpenFOAM. Stokes wave theories were applied to model the incoming waves in their work. As an extension to their study, the present work focuses on the investigation of the steep, non-breaking wave-induced soil response and liquefaction risk around gravity-based foundations of different geometries and incoming wave angles. The incoming waves in the present study are fully nonlinear and are modelled by fifth-order stream functions, according to Rienecker and Fenton (1981) and Fenton (1988). The nonlinear wave-structure interaction model is validated by comparing the free-surface elevation and the inline force with the experimental data conducted by Grue and Huseby (2002). The wave-seabed interaction model has been validated in the work of Li et al. (2018) by reproducing the laboratory experiment of Tsai and Lee (1995) of standing wave induced pore pressure under a vertical wall. A good agreement was obtained.

In the present work, gravity-based foundations with bottom slabs of cylindrical shape (a circular foundation) and hexagonal prismatic shape (a hexagonal foundation) are investigated. Fine-grained sandy seabed is modelled. At first, the soil consolidation states in the presence of the foundations in the static water are assessed. The initial effective stress in the soil is calculated. Then, the nonlinear waves are imposed. The nonlinear wave-induced pore pressure and liquefaction depth around the foundations are investigated for the hexagonal and circular foundations. For the hexagonal foundation, two different incoming wave angles with 0 degrees and 90 degrees, i.e. waves propagate to the hexagon corner and to the hexagon edge, are considered. The effect of the incoming wave angle on the liquefaction distribution around the hexagonal foundation is investigated.

2. Mathematical models and the coupling algorithm

2.1. Free surface wave model

The wave domain is governed by the incompressible Navier-Stokes equations including the continuity equation and the momentum equations.(1) $\mathrm{\nabla }\cdot \mathbf{u}=0$(1) (2) $\frac{\mathrm{\partial }\mathbf{u}}{\mathrm{\partial }t}+(\mathbf{u}\cdot \mathrm{\nabla })\mathbf{u}=-\frac{1}{\rho }\mathrm{\nabla }p+\mathbf{g}+\frac{1}{\rho }\mathrm{\nabla }\cdot \tau \tau$(2) where $\mathbf{u}$ is the velocity vector with three components in the x,y,and z directions respectively; $\mathbf{g}$ is the gravitational acceleration; ρ is the fluid density, p is the pressure and $\tau \tau$ (in Einstein notation ${\tau }_{ij}$) is the viscous stress tensor. For Newtonian fluid,(3) ${\tau }_{ij}=2\mu {\sigma }_{ij}$(3) where μ is the dynamic molecular viscosity ${\sigma }_{ij}$ is defined by(4) ${\sigma }_{ij}=\frac{1}{2}\left(\frac{\mathrm{\partial }{u}_i}{\mathrm{\partial }{x}_j}+\frac{\mathrm{\partial }{u}_j}{\mathrm{\partial }{x}_i}\right)$(4) where $i,j\in$ [1,2,3]. $u_i$ and $u_j$ are the velocity components in x, y and z direction respectively.

The equations are solved for the two immiscible fluids, i.e. air and water simultaneously. The fluids are tracked by using a scalar field α. $\alpha =0$ for air and $\alpha =1$ for water. The distribution of α is modelled by an advection equation(5) $\frac{\mathrm{\partial }\alpha }{\mathrm{\partial }t}+\mathrm{\nabla }\cdot \alpha \mathbf{u}+\mathrm{\nabla }\cdot \left[\alpha \right(1-\alpha \left){\mathbf{u}}_r\right]=0$(5) The last term on the left-hand side is a compression term, which limits the smearing of the interface, and ${\mathbf{u}}_r$ is a relative velocity which can be referred to Berberović et al. (2009).

Using α, the spatial variation in any fluid property can be expressed through the weighting(6) $\mathrm{\Phi }=\alpha {\mathrm{\Phi }}_{water}+(1-\alpha ){\mathrm{\Phi }}_{air}$(6) Φ is a fluid property, such as μ and ρ.

2.2. Linear elastic structure model

The structure domain is governed by a linear momentum balance equation,(7) $\mathrm{\nabla }\cdot \sigma \sigma =\mathrm{\nabla }\cdot [2\mu ϵϵ+\lambda \mathrm{t}\mathrm{r}(ϵ)ϵ)\mathbf{I}]=0$(7) and the isotropic linear elastic strain-displacement relations(8) $ϵϵ=1/2(\mathrm{\nabla }\mathbf{U}+\mathrm{\nabla }{\mathbf{U}}^T)$(8) where $\sigma \sigma$ is the stress tensor; $ϵϵ$ is the small strain tensor; $\mathbf{U}$ is the structural displacement vector consisting of three coupled components. $\mathrm{t}\mathrm{r}$ is the trace of an square matrix which is the sum of the elements on the main diagonal, $\mathrm{t}\mathrm{r}\left(ϵϵ\right)=\sum _{i=1}^n{ϵ}_{ii}$. μ and λ are the Lamé's coefficients of elastic material properties which relate to commonly used Young's modulus E and Poisson's ratio ν. The finite volume-based solution to the coupled equations can be referred to the work of Jasak and Weller (2000) and Li et al. (2018).

2.3. Anisotropic Biot's poro-elastic soil model

The seabed model is based on assumptions that the sandy seabed is horizontally flat, fully saturated, with anisotropic permeabilities in vertical and horizontal directions. The soil skeleton generally obeys Hooke's law with elastic properties, therefore the classical Biot's consolidation equations (Biot 1941) can be adopted. For anisotropic soil materials, the orthotropic elastic stress-strain relation can be expressed in a 6 × 6 matrix notation:(9) ${\sigma }^{\prime }{\sigma }^{\prime }=\left(\begin{array}{c}{\sigma }_{xx}^{\prime }\\ {\sigma }_{yy}^{\prime }\\ {\sigma }_{zz}^{\prime }\\ {\sigma }_{xy}\\ {\sigma }_{yz}\\ {\sigma }_{xz}\end{array}\right)=\left[\begin{array}{cccccc}{A}_{11}& A_{12}& A_{31}& 0& 0& 0\\ A_{12}& A_{22}& A_{23}& 0& 0& 0\\ A_{31}& A_{23}& A_{33}& 0& 0& 0\\ 0& 0& 0& A_{44}& 0& 0\\ 0& 0& 0& 0& A_{55}& 0\\ 0& 0& 0& 0& 0& A_{66}\end{array}\right]\left(\begin{array}{c}{ϵ}_{xx}\\ {ϵ}_{yy}\\ {ϵ}_{zz}\\ {ϵ}_{xy}\\ {ϵ}_{yz}\\ {ϵ}_{xz}\end{array}\right)=\mathbf{C}:ϵϵ$(9) where ${\sigma }^{\prime }{\sigma }^{\prime }$ is the effective stress tensor. In the present work, tension is taken as positive while compression is taken as negative. According to Demirdžić et al. (2000), the 9 independent coefficients $A_{ij}$ are calculated from Young's modulus $E_i$ and Poisson's ratio ${\nu }_{ij}$ and shear modulus $G_{ij}$. The calculation of the coefficients $A_{ij}$ can be referred to Li et al. (2018).

The soil domain is governed by a quasi-static momentum balance equation and a mass balance equation based on Darcy's law.

Quasi-static momentum balance equation for the soil mixture is expressed as:(10) $\mathrm{\nabla }\cdot \left[\mathbf{C}:\frac{1}{2}(\mathrm{\nabla }\mathbf{U}+(\mathrm{\nabla }\mathbf{U}{)}^T)\right]-\mathrm{\nabla }p=0$(10) where $\mathbf{U}$ is the soil (skeleton) displacement, p is the pore fluid pressure and $\mathbf{C}$ is the fourth-order elastic stiffness tensor.

Mass balance equation for the pore fluid is expressed as:(11) $\frac{n}{K^{\prime }}\frac{\mathrm{\partial }p}{\mathrm{\partial }t}-\frac{1}{{\gamma }_w}\mathrm{\nabla }\cdot (\mathbf{k}\cdot \mathrm{\nabla }p)+\frac{\mathrm{\partial }}{\mathrm{\partial }t}(\mathrm{\nabla }\cdot \mathbf{U})=0$(11) where n is the soil porosity, ${\gamma }_w$ is the specific weight of water in soil, and $\mathbf{k}$ is the diagonal permeability tensor with values of $k_x$, $k_y$ and $k_z$. The bulk modulus of the compressible pore flow $K^{\prime }$ is approximately computed by the formulation Vafai and Tien (1981):(12) $\frac{1}{K^{\prime }}=\frac{1}{K_w}+\frac{1-S_r}{p_a}$(12) where $S_r$ is the degree of soil saturation, $K_w$ is the bulk modulus of pure water and $p_a$ is the absolute pore water pressure.

The wave, structure and seabed domains are coupled by the data transfer at the common boundaries. A one-way coupling algorithm is adopted, as illustrated in Figure 1. The data is transferred in a single direction and the structure and seabed deformations do not alter the wave domain. This is because that the present study is focussed on the steep wave-induced momentary liquefaction, which means the phenomenon that investigated is an instantaneous phenomenon. Therefore, the time-dependent effect of wave-induced structure vibration on the wave and seabed domains and the effect of seabed deformation on the structure and wave domains are neglected.

Figure 1. One-way boundary coupling algorithm of wave-structure-seabed interaction.

3. Model validation

3.1. Validation of nonlinear wave-structure interaction

The present nonlinear wave-structure interaction model is validated against the experiment conducted by Grue and Huseby (2002). Their test of regular steep waves propagating towards a vertical circular cylinder (Figures 3(c,d) in Grue and Huseby (2002)) is reproduced. The free surface elevation and the inline force on the cylinder are compared between the present numerical simulation and the experimental data. The parameter setting for reproducing the experiment is shown in Table 1. D is the diameter (characteristic length) of the cylinder, $h_w$ is the water depth, ${\eta }_m$ is the wave crest elevation and k is the wave number. The fully nonlinear stream function wave theory is used in the numerical wave generation. A slip boundary condition is set at the cylinder surface and the seabed. The number of grid points for the free surface is 32 in the present study.

Table 1. Parameters of the steep wave propagation experiment of Grue and Huseby (2002).

D (m)	$h_w$ (m)	T (s)	$k{\eta }_m$	H (m)	${\eta }_m$
0.06	0.6	0.86	0.34	0.11	0.662

Figure 2 presents the comparison between the numerical results and the experimental measurements. Both the free surface elevation and the inline force are in good agreement (converging to 1–3% discrepancy as the waves develop). The experimental measurement contains second-order spurious waves from first-order wave generation. In the numerical simulation, the second-order spurious waves are avoided in the fully nonlinear stream function wave model, as pointed out by Paulsen et al. (2014).

Figure 2. Validation study of nonlinear wave propagation to a vertical structure by comparing to the experimental data of Grue and Huseby (2002).

3.2. Validation of wave-seabed interaction

The present wave-induced seabed response model has been validated in the work of Li et al. (2018). The experiment performed by Tsai and Lee (1995) of standing wave-induced soil response in the vicinity of a vertical wall was reproduced. It showed that the present numerical results are in good agreement with the measured data.

4. Model application

4.1. Numerical settings

The present work focuses on two types of gravity-based foundations. The geometries of the foundations are shown in Figure 3.

Figure 3. Geometries and meshes of the structures. (a) Foundation with cylindrical slab and (b) Foundation with hexagonal prismatic slab.

The entire multiphysics numerical system is established in the 3D Cartesian coordinate system. $z=0$ is located at the undisturbed free surface. x positive towards the wave propagation direction, y positive towards the side of the tank, and z positive upwards. The layout of the numerical system is shown in Figure 4. The wave inlet relaxation zone (Jacobsen et al. 2012) is set to be one wave length. The wave outlet relaxation zone is set to 1.25 wave lengths to absorb waves and ensure no reflection from the outlet boundary. The width of the wave tank (distance between the sides of the tank) is set to be two wave lengths.

Figure 4. A sketch of the numerical layout for the wave-structure-seabed interaction model.

Table 2 gives the parameter settings in the present work. The water depth is $h_w=10\mathrm{m}$, the wave height is $H=4\mathrm{m}$ and the wave period is $T=10\mathrm{s}$. According to Le Méhauté (2013), the present wave condition exceeds the modelling ranges of linear wave theory and Stokes wave theory. It should be modelled by using the fifth-order stream function theory (Dean 1965; Fenton 1988) for the steep non-breaking waves. The Keulegan-Carpenter number (KC number) (Journée and Massie 2000) for the simulated condition is 1.70, which is a small value that indicates the drag force which comes from the viscosity is negligible compared to inertial forces. Hence, turbulence modelling is not considered.

Table 2. Parameter settings for wave-structure-seabed interaction modelling in the present study.

Wave parameters
Wave height H (m)	4.0
Wave period T (s)	10.0
Water depth $h_w$ (m)	10.0
Wave length L (m)	98.2
KC number	1.70
Structure parameters
Characteristic length D (m)	19.0
Bottom slab height $h_b$ (m)	3.0
Shaft radius R (m)	2.25
Density ρ (kg/m${}^3$)	2400
Young's modules (N/m${}^2$)	$2.2\times {10}^{10}$
Poisson's ratio	0.2
Seabed parameters	Directional values
Young's modules (N/m${}^2$)	$E_x=1.2\times {10}^7$	$E_y=1.2\times {10}^7$	$E_z=2\times {10}^7$
Poisson's ratios	${\nu }_{xy}=0.2$	${\nu }_{yz}=0.24$	${\nu }_{zx}=0.4$
Shear modulus (N/m${}^2$)	$G_{xy}=5\times {10}^6$	$G_{yz}=1.2\times {10}^7$	$G_{zx}=1.2\times {10}^7$
Permeabilities (m/s)	$k_x=0.0005$	$k_y=0.0005$	$k_z=0.0001$
Saturation factor $S_r$	0.98
Porosity n	0.3

For the hexagonal foundations, two wave angles are considered in terms of waves come towards the hexagon corner (wave angle of 0 degrees) and waves come towards the hexagon edge (wave angle of 90 degrees), respectively. As shown in Figure 5, three points P1, P2 and P3 around the foundation are identified. The characteristic length of the foundation D is defined as the distance between P1 and P3.

Figure 5. Top views of waves incoming to the circular and hexagonal foundations. (a) circular foundation (wave angle $=0^{\circ }$). (b) hexagonal foundation (wave angle $=0^{\circ }$) and (c) hexagonal foundation (wave angle $={90}^{\circ }$).

4.2. Boundary conditions

Boundary conditions for the multiphysics model are specified as follows.

Wave domain:

At the inlet boundary, the velocity is specified as the input wave velocity. At the outlet boundary, the velocity is set to zero. At four sides of the numerical wave tank, the normal gradient of the pressure is set to zero, i.e. $\mathrm{\partial }p/\mathrm{\partial }n=0$. At the atmosphere, the pressure is set to be atmosphere pressure $p=p_0$. At the wave-structure interface and the wave-seabed interface, slip velocity boundary is set since the viscous effects at the walls are negligible.

Structure domain:

The structure transfers the dynamic wave loading to the seabed. At the wave-structure interface, the nonuniform and time-varying dynamic wave pressure is imposed on the structure. At the structure-seabed interface, the instantaneous structural stress is read by the soil solver to compute the compatible time-varying displacement gradient boundary for the seabed.

Seabed domain:

At the wave-seabed interface, the soil has zero effective shear stress while the pore pressure is equal to the dynamic wave pressure. At the structure-seabed interface, it is considered that the structure is impermeable. Thus, the pore pressure has zero normal gradient, according to Darcy's flow equation. At the bottom and lateral boundaries of the seabed, a slip boundary is applied for the soil skeleton. The normal pore pressure gradient is zero.

5. Results and discussions

5.1. Consolidation

Once the structure is placed on the seabed, the gravitational forces will induce the compression of the soil skeleton and the gradual dissipation of the excess pore pressure. The soil consolidation behaviour will further affect the potential risk of wave-induced soil liquefaction (Ulker et al. 2010; Sui et al. 2017; Zhao et al. 2017). Therefore, it is essential to do consolidation analysis before investigating the momentary liquefaction of the seabed around the foundations. In the present work, the consolidation process is simulated for each type of the structure. According to the 1D Terzaghi's consolidation theory, the time for completing 90% consolidation can be expressed as Wang (2017):(13) $t_{90}=T_v\frac{H_d^2}{c_v}$(13) where $H_d$ is the drainage distance of the layer, $T_v=0.848$ is the vertical consolidation time factor for 90% consolidation, $c_v$ is the consolidation coefficient calculated by(14) $c_v=\frac{2G{k}_z(1-\nu )}{{\gamma }_w(1-2\nu )}$(14) where $k_z$ is the vertical permeability and ${\gamma }_w$ is the bulk specific weight of the pore water. In the present study, the longest time for reaching the 90% consolidation state is estimated as 3230 s.

Figure 6 shows the vertical soil displacement $U_z$ after the completion of the consolidation in the presence of the circular foundation. A negative value of $U_z$ means the soil skeleton is compressed and moving downward. It is seen that around the foundations, the soil skeleton has been largely compressed. Compared to the circumstance without a structure, the presence of the foundation changes the surrounding soil displacement behaviour.

Figure 6. Initial vertical displacement (m) of the seabed after the completion of the consolidation with and without the circular foundation.

Figure 7 presents the distribution of the initial vertical effective stresses ${\sigma }_z^{\prime }$ in the soil after the consolidation process. A negative value of ${\sigma }_z^{\prime }$ indicates the compression of the soil skeleton. In the long-time consolidation process, the gravity forces gradually transfer to the supporting soil skeleton. It shows that the initial vertical effective stress below the foundation is amplified compared to the far field. Similar soil consolidation behaviour can be observed in the presence of the hexagonal foundation, as shown in Figure 8(a,b).

Figure 7. Initial vertical effective stress (N/m${}^2$) in the seabed after the completion of the consolidation with and without the circular foundation.

Figure 8. Initial vertical displacement (m) and initial vertical effective stress (N/m${}^2$) in the seabed after the completion of the consolidation with the hexagonal foundation. (a) Initial vertical displacement (m) of the seabed beneath the hexagonal foundation and (b) Initial vertical effective stress (N/m${}^2$) in the seabed beneath the hexagonal foundation.

A cross-sectional view of vertical soil effective stress at $z=-10.1\mathrm{m}$ (0.1 m below the seabed surface) is shown in Figure 9. It is observed that for the circular foundation, the initial vertical effective stress is evenly distributed in the surrounding soil. However, for the hexagonal foundation, the compressive stress at the corners is higher than that at the edges. This is because that the sharp corners of the hexagon cause stress concentration in the soil. The initial vertical effective stress ${\sigma }_z^{\prime }$ at the end of the consolidation process will be used as an initial condition to evaluate the momentary liquefaction risk.

Figure 9. A cross-section of the initial vertical effective stress (N/m${}^2$) of the seabed after the completion of the consolidation at 0.1 m below the seabed surface ($z=-10.1m$).

It is worthwhile to mention that the ‘pre-shear history’ of sand may have effect on its liquefaction resistance. Ye et al. (2015) studied the liquefaction risk of the sand with different initial shear stress ${\sigma }_{\theta z}$ (${\sigma }_{\theta z}$ correspond to clockwise shear load). They found that a high initial shear stress ratio (i.e. ${\sigma }_{\theta z}$ / mean effective stress) can speed up the liquefaction process. However, in the present study, the initial shear stress that is triggered by the consolidation process is negligible compared to the initial vertical effective stress. As shown in Figure 10, the magnitudes of initial shear stress ${\sigma }_{xz0}$ and initial vertical effective stress ${\sigma }_{z0}^{\prime }$ beneath the corner of the hexagonal foundation are compared. It is observed that the magnitude of ${\sigma }_{xz0}$ is negligible compared to ${\sigma }_{z0}^{\prime }$, which means that the initial shear stress ratio is very small. Therefore, the effect of ‘pre-shear history’ of sand on the liquefaction analysis is not considered in the present work.

Figure 10. Vertical distributions of initial shear stress and initial vertical effective stress after the consolidation process. It is observed that the magnitude of ${\sigma }_{xz0}$ is negligible compared to ${\sigma }_{z0}^{\prime }$. (a) Initial shear stress ${\sigma }_{xz0}$ beneath the corner of the hexagonal foundation and (b) Comparing the magnitudes of initial shear stress ${\sigma }_{xz0}$ and initial vertical effective stress ${\sigma }_{z0}^{\prime }$ beneath the corner of the hexagonal foundation.

5.2. Nonlinear wave-induced soil response

Nonlinear steep waves are imposed on the structure and seabed. Figure 11 presents the time histories of free surface elevation at the location A1 and the corresponding bottom wave pressure at B1. The locations of the points are shown in Figure 4. It is observed that the interaction between the waves and the foundations triggers strong nonlinearity of free surface elevation. It affects the bottom wave pressure, which shows a similar nonlinear effect.

Figure 11. Time histories of the free surface elevation and the corresponding bottom wave pressure. (a) circular foundation (wave angle $=0^{\circ }$). (b) hexagonal foundation (wave angle $=0^{\circ }$) and (c) hexagonal foundation (wave angle $={90}^{\circ }$).

The wave diffraction effect causes various pore pressure distributions at different locations around the foundation. According to Li et al. (2018), the soil beneath the foundation is shielded from the dynamic waves, so that the pore pressure underneath the foundation bottom has relatively small variations. Therefore, liquefaction is less likely to happen under the foundation bottom. The present work focuses on investigating the pore pressure distribution around the foundations. The vertical profiles of transient pore pressure at three typical locations surrounding the foundation are shown in Figure 12. P1 is at the upstream side of the foundation when waves first hit the surface of the foundations. P2 is at the location of 90 degrees from the upstream. P3 is at the downstream side of the foundation.

Figure 12. Transient pore pressure distribution at two time instants around the offshore foundations. (a) Transient pore pressure distribution around the circular foundation. (b) Transient pore pressure distribution around the hexagonal foundation with 0 degrees incoming waves and (c) Transient pore pressure distribution around the hexagonal foundation with 90 degrees incoming waves.

In Figure 12, vertical pore pressure profiles at time instants of $t/T=8$ and $t/T=8.5$ are presented. The vertical distance to the seabed surface $z^{\prime }$ is normalised by the seabed thickness $D_s$. Figure 12 shows that the magnitude of pore water pressure decreases rapidly from the seabed surface to approximately $z^{\prime }/D_s=-0.25$. Then, the magnitude remains almost unchanged until to the seabed bottom ($z=-D_s$). Comparing between the three curves in each sub-figure, it is observed that the pore pressure profiles at P1 (upstream) has the highest changes from the seabed surface to a depth at $z^{\prime }/D_s=-0.25$. This indicates that the pore pressure gradient at the upstream of the foundations are generally higher than that at the downstream of the foundations. Therefore, the liquefaction risk at the upstream can be relatively higher. The same phenomenon is found in the work of Lin et al. (2017). Their study showed that the liquefaction depth at the upstream side of the monopile is around 15–20% higher than that at the downstream side. However, this conclusion can only be drawn when the soil consolidation process is not considered. In fact, when the structure is placed on the seabed, the effective stress in the surrounding soil can significantly increase, which further influences the momentary liquefaction depth around the foundation (Ulker et al. 2010; Sui et al. 2017; Zhao et al. 2017).

Momentary liquefaction happens periodically under the wave trough. At this moment, upward seepage flows are generated. Figure 13 shows the vertical distribution of the maximum amplitude of negative pore pressure within a wave cycle. It is observed that, when incoming waves are perpendicular to the edge of the hexagonal foundation (Figure 13(c)), the vertical gradient of the pore pressure is the highest among the three cases. In Figure 13(b), incoming waves first hit the corner of the hexagon. The pore pressure and vertical pore pressure gradient at location P1 are smaller than that in Figure 13(c). This is because when waves hit the corner, the change of wave velocity in the wave propagating direction is smaller than when waves hit the edge. In the latter case, when waves hit the edge of the hexagon, the wave velocity component perpendicular to the edge becomes zero, so that the pressure is significantly increased. For the circular foundation, the change of the wave velocity is milder when waves diffract around the circular geometry, so that the pore pressure gradient around the circular foundation is generally smaller than that around the hexagonal foundations.

Figure 13. Vertical distribution of the amplitude of the negative pore pressure.

Based on the discussion above, the upstream of the hexagonal foundations can experience a higher pore pressure gradient than the circular foundation. If the initial effective stress caused by the seabed consolidation process is not considered, the liquefaction risk at the upstream of the hexagonal foundation should be higher than that at the upstream of the circular foundation.

5.3. Liquefaction

So far in the literature, various liquefaction criteria for investigating the momentary liquefaction have been proposed based on the either the effective stress (Okusa 1985; Tsai 1995), or the excess pore pressure (Zen and Yamazaki 1990; Jeng 1997). Ye (2012a) performed a comparison study of different liquefaction critera. For the liquefaction analysis considering the seabed consolidation around a marine structure, a modified criterion based on Zen and Yamazaki (1990) is recommended. The criterion is given in Equation (15(15) $p-p_b\ge {\sigma }_{z0}^{\prime }$(15) ).(15) $p-p_b\ge {\sigma }_{z0}^{\prime }$(15) where ${\sigma }_{z0}^{\prime }$ is the initial vertical effective stress induced by the gravitational forces from the consolidation process.

Maximum liquefaction depth in a wave cycle around the hexagonal foundations and the circular foundation are presented in Figure 14. It is seen that for the circular foundation, the liquefaction depth is approximately averagely distributed around the circular foundation. The upstream side has a sightly higher liquefaction depth than the downstream side. The reason is explained in Section 5.2.

Figure 14. Maximum liquefaction depth (m) in a wave cycle around the circular foundation (a) and the hexagonal foundations (b,c).

For the hexagonal foundations, the initial soil effective stress around the corners is higher than that around the edges, so that liquefaction is less likely to happen around the hexagon corners, as shown in Figure 14(b,c). However, when incoming waves propagate towards the hexagon corner (Figure 14(b)), the high pore pressure gradient can cause notable upward seepage flow. Therefore, liquefaction happens at the upstream corner in Figure 14(b). When incoming waves propagate ${90}^{\circ }$ to the hexagon edge, as shown in Figure 14(c), no liquefaction is observed at the hexagon corners.

In the real ocean environment, the waves are nonuniform and non-unidirectional, so that both the corners and edges of the hexagonal foundation can experience liquefaction and therefore need specific protections. Chang and Jeng (2014) studied the seabed protection methodology and suggested to replace the existing layers of the surrounding soil with materials of higher permeability to mitigate the liquefaction risk. It is noted that the process of soil liquefaction can be closely related to scour. In the experiment by Sumer et al. (2007), liquefaction and scour around a monopile were observed in one experiment. As discussed in Whitehouse (1998), since a liquefied bed has almost no shear strength, the sediment can be eroded more easily than a non-liquefied bed. In the engineering practice, gravity-based structures are usually fitted with a skirt to prevent the hydraulic process of scour channels penetrating underneath the structure.

6. Conclusions

The present study has investigated the soil consolidation, steep wave-induced soil response and the momentary liquefaction risk around gravity-based offshore foundations. Two different geometries are considered in terms of a circular foundation and a hexagonal foundation. For the hexagonal foundation, two different incoming wave angles are investigated in terms of waves come towards the hexagon corner and waves come towards the hexagon edge. A 3D FVM-based wave-structure-seabed interaction model has been applied. The nonlinear wave-structure interaction model and wave-seabed interaction model have been validated against existing experimental data. Good agreement has been obtained. An anisotropic poro-elastic model has been applied for the soil analysis. The following conclusions can be drawn from the present study:

1. The presence of the foundations on the seabed leads to the compression of the soil in the vicinity of the structure during the consolidation process. The initial vertical effective stress in the surrounding soil is increased. The slab geometry of the foundation affects the initial effective stress distribution. For the circular foundation, the initial vertical effective stress is evenly distributed in the surrounding soil. However, for the hexagonal foundation, the geometry triggers stress concentration in the soil around the corners. The compressive stress at the corners is much higher than that at the edges.

2. Generally, when waves propagate towards the foundations, the pore pressure gradient at the upstream of the foundations is higher than that at the downstream. Therefore, the liquefaction risk at the upstream of the foundation can be relatively higher. This is due to the wave diffraction effect.

3. For the circular foundation, the liquefaction depth is approximately averagely distributed around the circular foundation, except that the upstream side has a sightly higher liquefaction depth than the downstream side. For the hexagonal foundations, the distribution of the initial soil effective stress from the consolidation process affect the liquefaction zone distribution. The momentary liquefaction is less likely to happen around the hexagon corners due to high initial effective stress.

4. Wave direction affects the pore pressure and the liquefaction distribution around the hexagonal foundation. When incoming waves propagate towards the hexagon corner, the high pore pressure gradient can cause notable upward seepage flow. Although the initial effective stress at the hexagon corners is high, significant upward excess pore pressure can still trigger liquefaction around the upstream corner. However, when incoming waves propagate ${90}^{\circ }$ to the hexagon edge, no liquefaction is observed at the hexagon corners.

5. In the real ocean environment, the waves are nonuniform and non-unidirectional, so that both the corners and edges of the hexagonal foundation can experience liquefaction and therefore need liquefaction protections.

More experimental data is required before a final conclusion can be given. Meanwhile, the present method could be useful as an engineering tool for predicting the wave-induced soil response and momentary liquefaction risk around offshore foundations.

Acknowledgements

This study was supported in part with computational resources provided by the Norwegian Metacenter for Computational Science (NOTUR), under Project No: NN9372K.

Disclosure statement

No potential conflict of interest was reported by the authors.