Two-Dimensional Numerical Study of Seabed Response around a Buried Pipeline under Wave and Current Loading

The evaluation of the wave-induced seabed response around a buried pipeline has been widely studied. However, the analysis of seabed response around marine structures under the wave and current loadings are still limited. In this paper, an integrated numerical model is proposed to examine the wave and current-induced pore pressure generation, for instance, oscillatory and residual pore pressure, around a buried pipeline. The present wave–current model is based on the Reynolds-Averaged Navier–Stokes (RANS) equation with k-ε turbulence while Biot’s equation is adopted to govern the seabed model. Based on this numerical model, it is found that wave characteristics (i.e., wave period), current velocity and seabed characteristics such as soil permeability, relative density, and shear modulus have a significant effect on the generation of pore pressure around the buried pipeline.


Introduction
Submarine pipelines are frequently used in the transportation of hydrocarbons such as oil and natural gas from offshore platforms to onshore terminals and in the disposal of industrial and municipal waste in offshore engineering.Since these submarine pipelines play such an important role, it has become one of the significant concern in marine and geotechnical engineering.Generally, submarine pipelines are buried into the seabed or within a trench.In the ocean environment, as waves propagate over the ocean floor, they generate dynamic pressure fluctuations, which will further induce excess pore pressure and reduce effective stress within the seabed soil.As the excess pore pressure increases, the soil particles loses its strength and part of the seabed becomes unstable.In most of the cases, seabed instability is the primary cause that leads to marine structure failure.Therefore, it is essential for offshore engineers to evaluate the wave-induced seabed response.
To date, numerous studies on wave-induced seabed response with marine structures such as breakwaters [1], monopiles [2] and pipelines [3] or without marine structures [4,5] had already been carried out.Based on laboratory experiments and field studies [6], the wave-induced pore pressure classifies into two mechanisms, which are oscillatory and residual mechanisms, as shown in Figure 1.Oscillatory pore pressure generation is caused by the amplitude damping and phase lag in pore pressure whereas residual pore pressure generation is the build-up excess pore pressure, which caused by the contraction of soil resulted from cyclic loading [7].In the real ocean environment, water waves and current flows exist simultaneously.Hence, in addition to wave loadings, current loadings should also be taken into consideration.Recently, numerous researchers [8][9][10][11][12] have conducted studies to understand the interaction between wave, current, and seabed.In general, to understand the phenomenon of wave-current-seabed interaction (WCSI) and wave-current-seabed-structure interaction (WCSSI), three common approaches are widely adopted, which are the analytical solution, laboratory experiments, and numerical methods.For the study of WCSI, Zhang et al. [13] and Liao et al. [14] developed an analytical approximation to calculate the soil response of a porous seabed under the combination of wave and current loading.In many cases, to simplify the problems, seabed condition is assumed to be isotropic.In Zhou's [15] research, the soil response in anisotropic seabed with a buried marine pipeline was calculated using the analytical method.In the analytical approach, most researchers utilized the third approximation wave-current interactions as the governing equation to determine the dynamic pressure from the wave model, which eventually used as the boundary condition for the porous seabed analysis in the seabed model.
Besides analytical method, laboratory experiments have also been reported in the literature.Liu et al. [16] conducted a series of experimental studies using a one-dimensional cylinder.The purpose of a one-dimensional cylinder experiment is to determine the wave-induced oscillatory soil response and to obtain the vertical profile of the pore pressure distribution.Two or three-dimensional wave flume tank experiments [17,18] and centrifuge modeling [19,20] were conducted in the laboratory to study WCSSI phenomenon and the stability of a pipeline.There are advantages and limitations to one-dimensional cylinder experiments and two-dimensional wave flume tank experiments or centrifuge experiments.For a one-dimensional cylinder experiment, many points can be measured due to their thick soil layers, but only oscillatory pore pressure and the vertical pore pressure distribution is obtained, whereas for a two-dimensional experiment, the accumulation of build-up excess pore pressure can be analyzed, but fewer surface points are measured due to the shallow soil layers.Recently, Yang et al. [21] conducted a flume experiment using a new method of the grey value of water's image to study the initial movement of mud particles due to current loading.
As the problem becomes more complicated to be solved by analytical solution and laboratory experiments, numerical methods are usually performed.Several researchers [11] had conducted numerical simulations to study the wave-current induced seabed response without any presence of marine structures.With the inclusion of marine structures such as pipelines in the seabed, the soil responses vary under wave and current loading.Zhao et al. conducted a two-dimensional model to study the influences of pore pressure accumulations around the vicinity of a fully buried pipeline in the seabed [22] and protected in trench layer with partial backfills [23].Zhou et al. [24] simulated the pore pressures, effective stresses and liquefaction potential around a buried in a poroelastic seabed subjected to cnoidal wave loading.In both Zhao's and Zhou's studies, only wave loading takes into Figure 1.The conceptual sketch of the pore pressure mechanisms (Adapted from Jeng [7]).
In the real ocean environment, water waves and current flows exist simultaneously.Hence, in addition to wave loadings, current loadings should also be taken into consideration.Recently, numerous researchers [8][9][10][11][12] have conducted studies to understand the interaction between wave, current, and seabed.In general, to understand the phenomenon of wave-current-seabed interaction (WCSI) and wave-current-seabed-structure interaction (WCSSI), three common approaches are widely adopted, which are the analytical solution, laboratory experiments, and numerical methods.For the study of WCSI, Zhang et al. [13] and Liao et al. [14] developed an analytical approximation to calculate the soil response of a porous seabed under the combination of wave and current loading.In many cases, to simplify the problems, seabed condition is assumed to be isotropic.In Zhou's [15] research, the soil response in anisotropic seabed with a buried marine pipeline was calculated using the analytical method.In the analytical approach, most researchers utilized the third approximation wave-current interactions as the governing equation to determine the dynamic pressure from the wave model, which eventually used as the boundary condition for the porous seabed analysis in the seabed model.
Besides analytical method, laboratory experiments have also been reported in the literature.Liu et al. [16] conducted a series of experimental studies using a one-dimensional cylinder.The purpose of a one-dimensional cylinder experiment is to determine the wave-induced oscillatory soil response and to obtain the vertical profile of the pore pressure distribution.Two or three-dimensional wave flume tank experiments [17,18] and centrifuge modeling [19,20] were conducted in the laboratory to study WCSSI phenomenon and the stability of a pipeline.There are advantages and limitations to one-dimensional cylinder experiments and two-dimensional wave flume tank experiments or centrifuge experiments.For a one-dimensional cylinder experiment, many points can be measured due to their thick soil layers, but only oscillatory pore pressure and the vertical pore pressure distribution is obtained, whereas for a two-dimensional experiment, the accumulation of build-up excess pore pressure can be analyzed, but fewer surface points are measured due to the shallow soil layers.Recently, Yang et al. [21] conducted a flume experiment using a new method of the grey value of water's image to study the initial movement of mud particles due to current loading.
As the problem becomes more complicated to be solved by analytical solution and laboratory experiments, numerical methods are usually performed.Several researchers [11] had conducted numerical simulations to study the wave-current induced seabed response without any presence of marine structures.With the inclusion of marine structures such as pipelines in the seabed, the soil responses vary under wave and current loading.Zhao et al. conducted a two-dimensional model to study the influences of pore pressure accumulations around the vicinity of a fully buried pipeline in the seabed [22] and protected in trench layer with partial backfills [23].Zhou et al. [24] simulated the pore pressures, effective stresses and liquefaction potential around a buried in a poroelastic seabed subjected to cnoidal wave loading.In both Zhao's and Zhou's studies, only wave loading takes into consideration.Recently, Duan et al. performed a two-dimensional [25] and three-dimensional [26] numerical simulation to study the wave-current induced soil liquefaction around the buried pipeline.However, in their research, only the oscillatory soil response under a wave and current loading have been determined.The wave-current induced oscillatory soil response is influenced by the current characteristics (following or opposing currents), wave characteristics (i.e., wave height, wave period and water depth) and soil properties (i.e., backfill thickness, permeability, and degree of saturation).Relevant research on wave-current-seabed-structure interaction is still limited.Therefore, more studies still needed to be conducted to provide a better understanding of this phenomenon.
In the present study, the objective is to develop a numerical model to analyze the pore pressure accumulation around the vicinity of a fully buried pipeline when the seabed is subjected to wave and current loadings.This numerical analysis consists of two submodels, which are the wave-current model and the seabed model.In the wave-current model, RANS equation with k-ε turbulence governs the wave motion and current flow, whereas in the seabed model, Biot's quasi-static equations are employed to calculate soil response such as soil displacements, oscillatory pore pressure, residual pore pressure and the effective stress within the seabed.The influence of wave, current, and seabed characteristics such as current velocities, wave period, shear modulus, permeability, and relative density on the soil responses are presented and discussed in the later sections.

Methods
This numerical model is made up of two submodels: wave-current model and seabed model with a buried pipeline.These two submodels are integrated with one-way coupling method.For a one-way coupling method of fluid-seabed interactions [27], the wave pressure computed in the wave solver is introduced into the seabed solver as a boundary condition to solve for the seabed response.Figure 2 illustrates a wave-current-seabed-structure interaction (WCSSI), where x and z represents the Cartesian coordinates system, h is the seabed thickness, d is the water depth, l is the length of the computational domain, e represents the embedment depth of pipeline (i.e., the distance from the surface of the seabed to the center of the pipeline) and u is the initial current velocity.In this study, the length of the computational domain (l) is set to be three times the wavelength to ignore the influence of lateral boundaries.[25] and three-dimensional [26] numerical simulation to study the wave-current induced soil liquefaction around the buried pipeline.However, in their research, only the oscillatory soil response under a wave and current loading have been determined.The wave-current induced oscillatory soil response is influenced by the current characteristics (following or opposing currents), wave characteristics (i.e., wave height, wave period and water depth) and soil properties (i.e., backfill thickness, permeability, and degree of saturation).
Relevant research on wave-current-seabed-structure interaction is still limited.Therefore, more studies still needed to be conducted to provide a better understanding of this phenomenon.
In the present study, the objective is to develop a numerical model to analyze the pore pressure accumulation around the vicinity of a fully buried pipeline when the seabed is subjected to wave and current loadings.This numerical analysis consists of two submodels, which are the wave-current model and the seabed model.In the wave-current model, RANS equation with - turbulence governs the wave motion and current flow, whereas in the seabed model, Biot's quasi-static equations are employed to calculate soil response such as soil displacements, oscillatory pore pressure, residual pore pressure and the effective stress within the seabed.The influence of wave, current, and seabed characteristics such as current velocities, wave period, shear modulus, permeability, and relative density on the soil responses are presented and discussed in the later sections.

Methods
This numerical model is made up of two submodels: wave-current model and seabed model with a buried pipeline.These two submodels are integrated with one-way coupling method.For a one-way coupling method of fluid-seabed interactions [27], the wave pressure computed in the wave solver is introduced into the seabed solver as a boundary condition to solve for the seabed response.Figure 2 illustrates a wave-current-seabed-structure interaction (WCSSI), where x and z represents the Cartesian coordinates system, h is the seabed thickness, d is the water depth, l is the length of the computational domain, e represents the embedment depth of pipeline (i.e., the distance from the surface of the seabed to the center of the pipeline) and u is the initial current velocity.In this study, the length of the computational domain (l) is set to be three times the wavelength to ignore the influence of lateral boundaries.

Wave-Current Model
The wave and current interaction model is based on the Reynolds-Averaged Navier-Stokes equation to simulate the flow motion.For the incompressible fluids, the mass conservation equation and the momentum conservation equation can be expressed as below:

Wave-Current Model
The wave and current interaction model is based on the Reynolds-Averaged Navier-Stokes equation to simulate the flow motion.For the incompressible fluids, the mass conservation equation and the momentum conservation equation can be expressed as below: where x i represents the Cartesian coordinate, () i and () j represent the index tensor notion, u i is the ensemble mean velocity (m/s), p is the fluid pressure (Pa), ρ f is the fluid density (kg/m 3 ), t is the time, g is the acceleration (m/s 2 ), µ is the dynamic viscosity (Pa•s), and −ρ f u i u j is the Reynolds stress tensor.By applying the eddy-viscosity assumptions, Reynold stress term can be expressed as: where µ t is the turbulence viscosity (Pa•s), k is the turbulence kinetic energy (TKE, m 2 /s 2 ), and δ ij is the Kronecker delta.By substituting Equation (3) into Equation ( 2), the following equation can be obtained: where µ e f f = µ + µ t , which is the total effective viscosity (Pa•s).
For the prediction of a fully turbulent flow, the standard k-ε turbulence model [28] based on model transport equations for the turbulence kinetic energy and its rate of dissipation can be expressed as follows: where and C ε2 are the empirical coefficients determined from experiments.The empirical coefficients used in this study are based on previous studies [29]:

Seabed Model
As mentioned earlier, pore pressure generation classifies into two mechanisms: oscillatory pore pressure and residual pore pressure.The pore pressure generation can be expressed as follows: where u represents the wave and current-induced excess pore pressure at a specific point, u e is the oscillatory component whereas u p is the residual component that is further expressed as u p = 1 T T 0 udt where T denotes the wave period.

Oscillatory Soil Response
In this study, the seabed model is considered to be porous and hydraulically permeable.The soil skeleton and pore fluid are assumed to be compressible and obey the Hooke's Law.Therefore, Biot's theory [30] is utilized to govern the soil response.The mass conservation and force equilibrium equation can be written as follows: where represents the Laplace operator, γ w is the unit weight of water (N/m 3 ), n s is the soil porosity, k s is the soil permeability (m/s), G denotes the shear modulus of the seabed soil (N/m 2 ), ν is the Poisson's ratio, and u e represents the wave-current induced oscillatory soil response.The compressibility of the pore fluid (β s ) and the elastic volume strain of the soil matrix (ε s ) can be defined as: where K w is the true modulus of water (which is taken as 2 × 10 9 N/m 2 ), S r is the degree of saturation, and P wo is the absolute water pressure.u s and w s represent the soil displacement at xand z-direction respectively.

Residual Soil Response
For homogeneous and isotropic soil, the residual pore pressure can be derived from the onedimensional Biot's consolidation equation, which can be expressed as follows: where u p is the wave-current induced residual pore pressure, C v is the coefficient of consolidation and f represents the source term of the pore pressure generation.
With the presence of a buried pipeline in the seabed, the present study considers a two-dimensional plane strain problem.Hence, a slight modification is made to Equation ( 13), the governing equation for residual pore pressure [7,31] can be written as follows: where C v is the coefficient of consolidation and f (x, z, t) is the source term, a function of space and time.These two parameters can be defined as follows respectively: where u g is the generation of pore-water pressure, τ(x, z, t) is the shear-stress term, |τ(x,z,t)| α r σ 0 represents the induced cyclic shear-stress ratio, which defines the pore pressure accumulation, α r and β r are empirical coefficients that obtain from large-scale simple shear tests.Both the empirical coefficients have a correlational relationship with the relative density of soil (D r ), which can be expressed as follows [32]:

Boundary Condition
Equation ( 8), (9), and (10) represents the governing equations to determine the oscillatory pore pressure and soil displacements within the seabed and should be solved with the appropriate boundary conditions.Hence, boundary conditions need to be specified at the appropriate locations; at the surface of the seabed where wave interacts with the seabed, at the bottom of the seabed, at the lateral boundaries and the seabed-pipeline interface.
At the seabed surface, it is common that the vertical effective normal stress and shear stress vanishes, therefore the pore pressure at the seabed surface is assumed to be equal to wave pressure.σ sz = τ sxz = 0 and P s = P b at z = 0 (18) where P b denotes the wave pressure at the seabed surface, which can be obtained from the wave-current model.At the lateral boundaries of the seabed, it is said to be impermeable (i.e., zero flux) and zero horizontal displacements occur, i.e., Since the seabed is resting on a rigid and impermeable base, it is assumed that there are zero displacements and no vertical flow occurring at the bottom of the seabed.For infinite seabed thickness, The buried pipeline is assumed to be elastic and impermeable; hence, along the surface of the pipeline, the pore pressure to the normal gradient is assumed to be zero, which can be written as: It is also assumed that there is no relative displacement occurs between soil particles and pipeline at the interface,

Numerical Scheme
In this study, the integrated model consists of a wave-current model and a seabed model.The linear wave along with a steady current is simulated using computational fluid dynamics (CFD) software, FLOW-3D v11.2 (Flow Science, Inc., Santa Fe, New Mexico, USA).From the wave-current simulation, wave pressure acts upon on the seabed surface can be obtained.Besides wave pressure, in FLOW-3D, the Volume of Fluid (VOF) method [33] is employed to get the free surface wave elevation.Then, the wave pressure on the seabed surface is introduced into the seabed model as a boundary condition.
As mentioned in Paulsen et al.'s report [34], an increase in spatial resolution leads to a reduction in errors.Therefore, the spatial resolution of an incident wave is recommended to be 15 points per wave height (p.p.w.h.).In this research, a total of 250,000 real cells with an element size of 0.2 m in both xand z-direction (i.e., 20 points per wave height) are constructed for the simulation of the wave-current model.
The seabed model is the simulation of the seabed and buried pipeline using a finite element analysis software, COMSOL Multiphysics v5.3a (COMSOL Inc., Burlington, MA, USA).COMSOL Multiphysics is a finite-element-method (FEM) software utilized to model and solve various types of scientific and engineering problems [35].Figure 3 shows a comparison between three different levels of predefined mesh range from COMSOL Multiphysics; i.e., fine (0.09 m-15.9 m), finer (0.0375 m-11.1 m) and extra fine (0.0225 m-6 m) mesh size.It can be observed in Figure 3a that smaller mesh sizes lead to an increase in the number of elements and longer computational time in solving the simulation.However, there is no visible difference in the result as shown in Figure 3b.The residual pore pressure over a specified period as shown in Figure 3b is taken from a point at x = 50 m and z = −5 m of the seabed.Hence, in this research, the seabed submodel consists of 2367 triangular elements with predefined finer mesh, which has mesh sizes range from 0.0375 m (domain near the pipeline) to 11.1 m (domain further away from the pipeline).
Both softwares, i.e.FLOW-3D and COMSOL Multiphysics, are running on Intel®Core (TN) 17-6700 CPU @ 3.40GHz with available memory of 8 GB.After building and simulating the model in COMSOL Multiphysics, MATLAB code is used to perform the loop simulation for the N -th number of computational time we intend to study.The MATLAB version used in this research is R2014a (The MathWorks, Inc., Natick, Massachusetts, USA).The MATLAB script is created to perform a series of iterative calculation with transient studies.The solution from the previous transient study is then set as the initial condition for the next transient study.For instance, at the first step during the time interval Δt = 1/nT, the wave load is introduced into the seabed submodel to calculate the soil response and study 1 is created.Then the solution from study 1 is set as the initial condition for study 2. The following steps repeat until the desired computational time is reached.The number of loop simulations corresponds to the number of wave load extracted from the wave-current submodel.
In this study, a total of 200 loop simulations were performed for a duration of 200 s, i.e., time interval of 1 s, to generate the dynamic soil responses for a specific model.The computational time is chosen as 200 s because at t = 200 s, the soil around the buried pipeline has begun to liquefy even though the oscillatory pore pressure and residual pore pressure from the simulations have not reached its steady state.The liquefaction criterion proposed Zen and Yamazaki [6] adopted to evaluate the liquefaction potential is expressed as: where γs is the unit weight of seabed soil, γw is the unit weight of water, z represents depth of a specific point in the seabed, Ps is the wave-current-induced pore pressure, and Pb is the wave pressure on the seabed surface.The parameters utilized in the simulation of the numerical model are listed in Table 1.After building and simulating the model in COMSOL Multiphysics, MATLAB code is used to perform the loop simulation for the N -th number of computational time we intend to study.The MATLAB version used in this research is R2014a (The MathWorks, Inc., Natick, MA, USA).The MATLAB script is created to perform a series of iterative calculation with transient studies.The solution from the previous transient study is then set as the initial condition for the next transient study.
For instance, at the first step during the time interval ∆t = 1/nT, the wave load is introduced into the seabed submodel to calculate the soil response and study 1 is created.Then the solution from study 1 is set as the initial condition for study 2. The following steps repeat until the desired computational time is reached.The number of loop simulations corresponds to the number of wave load extracted from the wave-current submodel.
In this study, a total of 200 loop simulations were performed for a duration of 200 s, i.e., time interval of 1 s, to generate the dynamic soil responses for a specific model.The computational time is chosen as 200 s because at t = 200 s, the soil around the buried pipeline has begun to liquefy even though the oscillatory pore pressure and residual pore pressure from the simulations have not reached its steady state.The liquefaction criterion proposed Zen and Yamazaki [6] adopted to evaluate the liquefaction potential is expressed as: where γ s is the unit weight of seabed soil, γ w is the unit weight of water, z represents depth of a specific point in the seabed, P s is the wave-current-induced pore pressure, and P b is the wave pressure on the seabed surface.The parameters utilized in the simulation of the numerical model are listed in Table 1.

Wave Characteristics
In the wave-current model simulation, the wave pressure and free surface elevation are simulated based on the wave characteristics that are input into the solver.The wave characteristics such as wave height (H), wave period (T), water depth (d) and the presence of current (u) influence the wave outcomes.The influence of the current velocities will be further discussed in the following section.

Influence of Current Velocities
Figure 4 shows the influence of different current flow on wave at a specific location over time.The information is extracted from the midpoint of the seabed surface within the computational domain (i.e., x = 150 m).Seabed surface midpoint (x = 150 m) is chosen as the specific point to be analyzed because wave trains are assumed to be more stable as there might be a minor disturbance at the start and end point of the computational domain (i.e., x = 0 m and x = 300 m respectively).At the starting point of the domain (i.e., x = 0 m), wave and current flow are first introduced.A wave absorber is placed after the computational domain (i.e., x = 300 m), which is for the dissipation of wave energy.From Figure 4, it is observed that as the current flow increases, wave pressure also increases.

Wave Characteristics
In the wave-current model simulation, the wave pressure and free surface elevation are simulated based on the wave characteristics that are input into the solver.The wave characteristics such as wave height (H), wave period (T), water depth (d) and the presence of current (u) influence the wave outcomes.The influence of the current velocities will be further discussed in the following section.

Influence of Current Velocities
Figure 4 shows the influence of different current flow on wave at a specific location over time.The information is extracted from the midpoint of the seabed surface within the computational domain (i.e., x = 150 m).Seabed surface midpoint (x = 150 m) is chosen as the specific point to be analyzed because wave trains are assumed to be more stable as there might be a minor disturbance at the start and end point of the computational domain (i.e., x = 0 m and x = 300 m respectively).At the starting point of the domain (i.e., x = 0 m), wave and current flow are first introduced.A wave absorber is placed after the computational domain (i.e., x = 300 m), which is for the dissipation of wave energy.From Figure 4, it is observed that as the current flow increases, wave pressure also increases.

Seabed Characteristics
In COMSOL Multiphysics, several simulations were conducted to simulate the wave-current induced soil response, i.e., oscillatory pore pressure and residual pore pressure, and the parameters used are as listed in Table 1, and the current velocity is taken as 0.5m/s in all the following simulations in Section 3.2.The generation of pore pressure in the seabed are sensitive to the seabed characteristics, for instance, soil permeability ( ), shear modulus of the seabed (G), relative density ( ) and degree of saturation ( ).The degree of sensitivity of these parameters on the pore pressure generation is discussed in the following section.

Effects of Soil Permeability
Soil permeability or also known as the hydraulic conductivity is a soil property, which allows the seepage of fluids to pass through its interconnected void spaces.Depending on the soil types, soil permeability can vary on a range from 1.0×10 −12 to 1.0×10 −2 m/s.In this analysis, three different values of soil permeability (i.e., k = 1.0×10 −2 , 1.0×10 −3 and 1.0×10 −4 m/s) are evaluated.Figure 5

Seabed Characteristics
In COMSOL Multiphysics, several simulations were conducted to simulate the wave-current induced soil response, i.e., oscillatory pore pressure and residual pore pressure, and the parameters used are as listed in Table 1, and the current velocity is taken as 0.5m/s in all the following simulations in Section 3.2.The generation of pore pressure in the seabed are sensitive to the seabed characteristics, for instance, soil permeability (k s ), shear modulus of the seabed (G), relative density (D r ) and degree of saturation (S r ).The degree of sensitivity of these parameters on the pore pressure generation is discussed in the following section.

Effects of Soil Permeability
Soil permeability or also known as the hydraulic conductivity is a soil property, which allows the seepage of fluids to pass through its interconnected void spaces.Depending on the soil types, soil permeability can vary on a range from 1.0 × 10 −12 to 1.0 × 10 −2 m/s.In this analysis, three different values of soil permeability (i.e., k = 1.0 × 10 −2 , 1.0 × 10 −3 and 1.0 × 10 −4 m/s) are evaluated.Figure 5a,b illustrate the distribution of wave and current induced oscillatory pore pressure and residual pore pressure at the point x = 50 m and z = −10 m respectively under different soil permeability conditions.The results expressed in Figure 5 are for the case in which shear modulus (G) and relative density (D r ) set as 5.0 × 10 6 N/m 2 and 0.5, respectively.
From Figure 5a, it can be observed that changes in the soil permeability have minimal influence on the value for oscillatory pore pressure.However, as noted in Figure 5b, seabed with low permeability (i.e., k s = 1.0 × 10 −4 m/s) tends to generate a higher value of residual pore pressure.It is because water cannot dissipate efficiently from low permeable soils and eventually resulted in the build-up of excess pore pressure.The results expressed in Figure 5 are for the case in which shear modulus (G) and relative density ( ) set as 5.0×10 6 N/m 2 and 0.5, respectively.From Figure 5(a), it can be observed that changes in the soil permeability have minimal influence on the value for oscillatory pore pressure.However, as noted in Figure 5(b), seabed with low permeability (i.e., ks = 1.0×10 −4 m/s) tends to generate a higher value of residual pore pressure.It is because water cannot dissipate efficiently from low permeable soils and eventually resulted in the build-up of excess pore pressure.

Effects of Shear Modulus
Figure 6 illustrates the variations of wave-current induced pore pressure at x = 50 m and z = −10 m under two different values of shear modulus (i.e., G = 5.0×10 6 and 1.5×10 7 N/m 2 ).In this section, the soil permeability and relative density are fixed at 0.0001 m/s and 0.5 respectively.Shear modulus is one of the soil properties that used to describe the tendency of an object deforms in shape at constant volume when acted upon by the opposing forces.Consider that the soil to be elastic, shear modulus (G) has a relationship with Young's modulus (E) and Poisson's ratio (v), which can be calculated from: Therefore, with a fixed value of Poisson's ratio (v) at 0.35, Young's modulus can be obtained from Equation (24) as 13.5 MPa and 40.5 MPa for shear modulus (G) of 5.0×10 6 N/m 2 and 1.0×10 7 N/m 2 respectively.The soil becomes denser when shear modulus and Young's modulus increases.From the reference table below (Table 2), soil with Young's modulus of 13.5 MPa is said to be loose sand while soil with Young's modulus of 40.5 MPa is said to be medium dense sand.When loose, saturated

Effects of Shear Modulus
Figure 6 illustrates the variations of wave-current induced pore pressure at x = 50 m and z = −10 m under two different values of shear modulus (i.e., G = 5.0 × 10 6 and 1.5 × 10 7 N/m 2 ).In this section, the soil permeability and relative density are fixed at 0.0001 m/s and 0.5 respectively.Shear modulus is one of the soil properties that used to describe the tendency of an object deforms in shape at constant volume when acted upon by the opposing forces.Consider that the soil to be elastic, shear modulus (G) has a relationship with Young's modulus (E) and Poisson's ratio (v), which can be calculated from: Therefore, with a fixed value of Poisson's ratio (v) at 0.35, Young's modulus can be obtained from Equation (24) as 13.5 MPa and 40.5 MPa for shear modulus (G) of 5.0 × 10 6 N/m 2 and 1.0 × 10 7 N/m 2 respectively.The soil becomes denser when shear modulus and Young's modulus increases.From the reference table below (Table 2), soil with Young's modulus of 13.5 MPa is said to be loose sand while soil with Young's modulus of 40.5 MPa is said to be medium dense sand.When loose, saturated soil is subjected to shear force, the soil particles tend to rearrange themselves into a denser manner,

Effects of Relative Density
Relative density (D r ) is a soil parameter that commonly used to indicate the in-situ denseness or looseness of granular soil [37].It is defined as the ratio of the difference between the void ratios of a cohesionless soil in its loosest state and existing natural state to the difference between its void ratio in the loosest and densest states, which can be formulated as follows: D r = e max − e e max − e min (25) As seen in Equation ( 17), when there is a change in the relative density of the soil, both the empirical coefficient α r and β r in the source term f (x, z, t) varies, which eventually leads to different values of residual pore pressure.Figure 6 shows the generation of oscillatory and residual pore pressure due to the influence of various relative density (D r = 0.2, 0.3 and 0.5) over a certain period at point x = 50 m and z = -10 m.In this section, the soil permeability (k) and shear modulus (G) are set as 1.0 × 10 −4 m/s and 5.0 × 10 6 N/m 2 respectively while other parameters remain unchanged as shown in Table 1. Figure 7a shows that oscillatory pore pressure is not affected by the change in relative density.However, the generation of residual pore pressure varies drastically as shown in Figure 7b.Higher accumulation of pore pressure is generated when the relative density is of a smaller value, i.e., looser soil.This concludes that loose sand tends to generate higher residual pore pressure (same explanation as Section 3.2.1).Relative density ( ) is a soil parameter that commonly used to indicate the in-situ denseness or looseness of granular soil [37].It is defined as the ratio of the difference between the void ratios of a cohesionless soil in its loosest state and existing natural state to the difference between its void ratio in the loosest and densest states, which can be formulated as follows: As seen in Equation ( 17), when there is a change in the relative density of the soil, both the empirical coefficient  and  in the source term (, , ) varies, which eventually leads to different values of residual pore pressure.Figure 6 shows the generation of oscillatory and residual pore pressure due to the influence of various relative density ( = 0.2, 0.3 and 0.5) over a certain period at point x = 50 m and z = -10 m.In this section, the soil permeability (k) and shear modulus (G) are set as 1.0×10 -4 m/s and 5.0×10 6 N/m 2 respectively while other parameters remain unchanged as shown in Table 1. Figure 7(a) shows that oscillatory pore pressure is not affected by the change in relative density.However, the generation of residual pore pressure varies drastically as shown in Figure 7(b).Higher accumulation of pore pressure is generated when the relative density is of a smaller value, i.e., looser soil.This concludes that loose sand tends to generate higher residual pore pressure (same explanation as Section 3.2.1).

Around the Vicinity of the Pipeline
In this section, we will study the dynamic soil response around the vicinity of the buried pipeline when the seabed is subjected to wave and current loading.The parameters utilized in the comparison of soil responses such as displacement and pore pressure around the buried pipeline under different current velocity are stated in Table 1, however, the shear modulus, permeability and relative density

Around the Vicinity of the Pipeline
In this section, we will study the dynamic soil response around the vicinity of the buried pipeline when the seabed is subjected to wave and current loading.The parameters utilized in the comparison of soil responses such as displacement and pore pressure around the buried pipeline under different current velocity are stated in Table 1, however, the shear modulus, permeability and relative density of the seabed in the simulation are set at a fixed value of 5.0 × 10 6 N/m 2 , 1.0 × 10 -4 m/s and 0.2 respectively.
As shown in the sketch of wave-current-seabed-pipeline interaction (Figure 2), the pipeline is located at x = 150 m and is embedded 3 m into the ground (Embedment depth is measured from seabed surface to the center of the pipeline).Hence, the analysis is taken at the point x = 150 m at t = 200 s.The graphs in Figure 8 show a comparison between various current velocity ranges from 0 m/s to 0.5 m/s in a gradient of 0.25, which are illustrated in terms of different soil response such as oscillatory pore pressure, residual pore pressure and displacement.The negative values of oscillatory pore pressure in Figure 8a is due to the wave trough phase at that specific time.It can be observed that with lower current values, the oscillatory pore pressure is higher.Figure 8b illustrates the residual pore pressure near the buried pipeline.It can be seen that at the surface of the seabed, residual pore pressure equals to zero.As it goes into the seabed closer to the top of the pipeline, the values increase gradually.The residual pore pressure tends to be higher at the bottom of the pipeline than the top of the pipeline.Figure 8c presents the displacement of seabed when wave and current loading apply on the seabed.It is observed that displacement is higher at the mid-depth of the seabed and gradually decreases to zero as it gets deeper into the seabed.It indicates that at a greater depth, the wave and current loading does not affect the movement of the soil particles.
of the seabed in the simulation are set at a fixed value of 5.0×10 6 N/m 2 , 1.0×10 -4 m/s and 0.2 respectively.
As shown in the sketch of wave-current-seabed-pipeline interaction (Figure 2), the pipeline is located at x = 150 m and is embedded 3 m into the ground (Embedment depth is measured from seabed surface to the center of the pipeline).Hence, the analysis is taken at the point x = 150 m at t = 200 s.The graphs in Figure 8 show a comparison between various current velocity ranges from 0 m/s to 0.5 m/s in a gradient of 0.25, which are illustrated in terms of different soil response such as oscillatory pore pressure, residual pore pressure and displacement.The negative values of oscillatory pore pressure in Figure 8(a) is due to the wave trough phase at that specific time.It can be observed that with lower current values, the oscillatory pore pressure is higher.Figure 8(b) illustrates the residual pore pressure near the buried pipeline.It can be seen that at the surface of the seabed, residual pore pressure equals to zero.As it goes into the seabed closer to the top of the pipeline, the values increase gradually.The residual pore pressure tends to be higher at the bottom of the pipeline than the top of the pipeline.Figure 8(c) presents the displacement of seabed when wave and current loading apply on the seabed.It is observed that displacement is higher at the mid-depth of the seabed and gradually decreases to zero as it gets deeper into the seabed.It indicates that at a greater depth, the wave and current loading does not affect the movement of the soil particles.Figure 9 displays the distribution of wave and current induced seabed responses; oscillatory pore pressure, residual pore pressure, and soil displacements,  and  around the periphery of the buried pipeline at a duration of 200 s (t = 200 s). Figure 9(a) and 9(b) show the horizontal and vertical soil displacement around the embedded pipeline respectively.There is no significant Figure 9 displays the distribution of wave and current induced seabed responses; oscillatory pore pressure, residual pore pressure, and soil displacements, u s and w s around the periphery of the buried pipeline at a duration of 200 s (t = 200 s). Figure 9a,b show the horizontal and vertical soil displacement around the embedded pipeline respectively.There is no significant difference in the soil displacement at the top and bottom of the pipeline.At t = 200 s, a wave trough reaches the proximity of the buried pipeline; therefore, a negative value of oscillatory pore pressure is observed as shown in Figure 9c.As shown in Figure 9d, the residual pore pressure increases gradually with depth because pore water pressure can dissipate efficiently at the seabed surface.Meanwhile, as it goes deeper into the seabed and pipeline, excess pore pressure accumulates which results in an increase in residual pore pressure at the seabed bottom.
difference in the soil displacement at the top and bottom of the pipeline.At t = 200 s, a wave trough reaches the proximity of the buried pipeline; therefore, a negative value of oscillatory pore pressure is observed as shown in Figure 9(c).As shown in Figure 9(d), the residual pore pressure increases gradually with depth because pore water pressure can dissipate efficiently at the seabed surface.Meanwhile, as it goes deeper into the seabed and pipeline, excess pore pressure accumulates which results in an increase in residual pore pressure at the seabed bottom.

Conclusion
This research proposes a two-dimensional model with a turbulence closure scheme (  - turbulence) to investigate the dynamic seabed responses around the vicinity of a buried pipeline under combined wave and current loading.However, according to Alberello et al. [38], there has no visible effect on the two-dimensional model with or without turbulence closure scheme whereas a substantial difference can be seen in a three-dimensional with turbulence model.Therefore, the current two-dimensional model is expected to expand into a three-dimensional model with a turbulence closure scheme in the future.
The effects of current, wave and soil characteristics on the wave-current-induced soil responses are examined.According to the numerical results presented above, the following conclusions can be drawn: (1) In the analysis of wave-seabed-structure interaction, the current should be taken into consideration-an increase in the current flow results in an increase in wave pressure.When

Conclusions
This research proposes a two-dimensional model with a turbulence closure scheme (k-ε turbulence) to investigate the dynamic seabed responses around the vicinity of a buried pipeline under combined wave and current loading.However, according to Alberello et al. [38], there has no visible effect on the two-dimensional model with or without turbulence closure scheme whereas a substantial difference can be seen in a three-dimensional with turbulence model.Therefore, the current two-dimensional model is expected to expand into a three-dimensional model with a turbulence closure scheme in the future.
The effects of current, wave and soil characteristics on the wave-current-induced soil responses are examined.According to the numerical results presented above, the following conclusions can be drawn: (1) In the analysis of wave-seabed-structure interaction, the current should be taken into consideration-an increase in the current flow results in an increase in wave pressure.When wave pressure increases, the oscillatory pore pressure tends to increase with increasing current velocity.(2) Soil permeability governs the seepage of fluid passing through or flowing out of the seabed.Low permeability, i.e., pore fluids cannot dissipate efficiently, resulted in higher residual pore pressure due to the increase in the buildup of excess pore pressure.(3) Shear modulus has a relationship with Young's modulus and Poisson's ratio, which describe the rigidity of the seabed.Keeping the Poisson's ratio as a constant value, a higher value of shear modulus generates a higher value of Young's modulus, which represents a denser soil.As presented above, loose sand tends to produce a higher value of residual pore pressure.(4) Relative density controls the empirical coefficients α r and β r in source term, which affects the generation of residual pore pressure.It concludes that a smaller value of relative density results in a higher value of residual pore pressure.However, there is no visible difference in the oscillatory pore pressure.

JFigure 3 .
Figure 3.Comparison between different mesh sizes for COMSOL Multiphysics in terms of (a) the number of elements and the computational time, and (b) residual pore pressure,  at a duration of 100 s.

Figure 3 .
Figure 3.Comparison between different mesh sizes for COMSOL Multiphysics in terms of (a) the number of elements and the computational time, and (b) residual pore pressure, u p at a duration of 100 s.

Figure 3 .
Figure 3. Variation of wave pressure,  at the seabed bottom for different current velocities.
(a) and 5(b) illustrate the distribution of wave and current induced oscillatory pore pressure and residual pore pressure at the point x = 50 m and z = −10 m respectively under different soil permeability conditions.

Figure 4 .
Figure 4. Variation of wave pressure, p at the seabed bottom for different current velocities.

Figure 4 .
Figure 4. Variations of (a) oscillatory pore pressure,  and (b) residual pore pressure,  at a specific point for different soil permeability.

Figure 5 .
Figure 5. Variations of (a) oscillatory pore pressure, u e and (b) residual pore pressure, u p at a specific point for different soil permeability.

Figure 5 .
Figure 5. Variations of (a) oscillatory pore pressure,  and (b) residual pore pressure,  at a specific point under the different shear modulus.

Figure 6 .
Figure 6.Variations of (a) oscillatory pore pressure, u e and (b) residual pore pressure, u p at a specific point under the different shear modulus.

Figure 6 .
Figure 6.Variation of (a) oscillatory pore pressure,  and (b) residual pore pressure,  _p at a specific point under the different relative density.

Figure 7 .
Figure 7. Variation of (a) oscillatory pore pressure, u e and (b) residual pore pressure, u p _p at a specific point under the different relative density.

Figure 8 .
Figure 8. Variation in (a) oscillatory pore pressure, u e , (b) residual pore pressure, u p and (c) displacement, u s under various current velocity.

Figure 8 .
Figure 8. Distribution of seabed responses such as (a) horizontal displacement  , (b) vertical displacement  , (c) oscillatory pore pressure,  , and residual pore pressure  around the vicinity of the buried pipeline due to wave and current loading when current velocity is set at  = 0.5 m/s.

Figure 9 .
Figure 9. Distribution of seabed responses such as (a) horizontal displacement u s , (b) vertical displacement w s , (c) oscillatory pore pressure, u e , and residual pore pressure u p around the vicinity of the buried pipeline due to wave and current loading when current velocity is set at v c = 0.5 m/s.

Table 1 .
Input data for the numerical simulation.
Commented [M1]: Please add the head

Table 1 .
Input data for the numerical simulation.