#### Abstract

Monitored breakwater settlements taken from an actual breakwater structure over an extended period of time (more than five years) were analyzed. The analysis revealed that the waves clearly affect the settlement of the breakwater, especially during high wave conditions such as typhoons. Breakwater settlement is caused by a decrease of effective stress of seabed during partial liquefaction due to wave-induced cyclic loads, which occurs due to an increase in excess pore pressure and the combination of oscillatory and residual pore water pressures. A new combined numerical model was suggested that allows the storm wave-induced seabed settlement underneath the caisson breakwater to be examined qualitatively. The technique uses a combined wave model (2D-NIT) and soil model (FLIP). The dynamic wave load calculated by the 2D-NIT was used as the input data for the soil model. This soil model can simulate both oscillatory and residual pore water pressures at the same time. There is a different feature to other previous studies adopting similar techniques.

#### 1. Introduction

Two key calculations are required to design a breakwater on sand layers: an assessment of the ultimate bearing capacity and an estimation of the settlement under working loads. Settlements are considered tolerable if they do not impair the functionality or serviceability of the foundations or the supported superstructures under the design loads [1–3].

Figure 1 shows the schematic diagram from the settlement data measured at the edge of the caisson breakwater placed at the west breakwater in Jeju’s outer harbor. As shown in this figure, a settlement occurs not only immediately after loading but also continually over time. Soil strata from the top layer at this site contain layers of marine sediment, sandy soil, and weathered rock. According to the soil profile, there is not a soft clay layer at this site. Kim et al. [4] studied this long-term compressional settlement problem. They conducted compression tests to assess the compression characteristics of the sands and found that the settlement issues are due to the high compressibility and particle crushing and shattering of the sand materials at the site, about 50% of which includes carbonate sands derived from shell. The carbonate sand is weaker than silicate sand and has many pores on the particle’s surface and an angular shape [5].

**(a)**

**(b)**

Another interesting phenomenon in Figure 1 is the rapid and abruptly large settlement at intervals *A* and *B*. This settlement is not related to the carbonate sands because its magnitude is fairly large, and it occurs over a very short period of time. The settlement at interval *A* occurred after Typhoon Maemi passed, as shown in Figure 1. Typhoon Maemi had a maximum wave height of 11.38 m and a maximum wave period of 7.8 sec. It was the first big storm after the caisson was placed on the site. Interval *B* indicates another big wave that attacked the caisson, similar to Typhoon Maemi. Some settlements also occurred here, but the magnitude was smaller than that of interval *A*. Based on the settlement record data, the structure built on sand was also exposed to wave-induced settlements. This phenomenon also has been captured by previous research done by Ye et al. [6].

Many numerical and experimental studies for wave-induced soil responses around breakwater have been carried out since the 1970s. Zen and Yamazaki [7], Mase et al. [8], Mizutani et al. [9], Ulker et al. [10], Li and Jeng [11], and Jeng et al. [12] used a poroelastic model for the seabed, and the dynamic Biot equation [13], known as “u-p” approximation, was proposed by Zienkiewicz et al. [14] to govern the dynamic response of the porous medium under wave loading. Ye and Jeng had conducted milestone work in this field. Ye et al. [15] developed a semicoupled numerical mode FSSI-CAS 2D to investigate the complicated interaction between ocean waves, offshore structures, and their loose seabed foundations. They further extended the 2D package to form FSSI-CAS 3D for the same purpose [16].

Oumeraci [17] analyzed many breakwater failures and concluded that the partial or complete liquefaction risk must be considered in the case of wave loading, which is especially relevant with a vertical breakwater or an offshore gravity platform. de Groot et al. [18,19] indicated that the foundation of the vertical breakwater, or an offshore gravity platform built in sandy soil, may fail due to wave-induced partial liquefaction. Kudella [20] conducted large-scale model experiments in a wave flume to study the generation of transient, instantaneous, and residual pore pressure in a seabed beneath a caisson breakwater that has been subjected to wave loads. These experiments showed that the residual pore pressure is essentially generated by the caisson’s movements when breaking the wave loads and that they are closely related to residual soil deformation. These previous studies clearly show the effect of waves on seabed deformation underneath the breakwater.

Like previous studies done by other researchers, this study focused on the dynamic response of a sandy seabed and the liquefaction in the sand seabed’s foundation under typhoon-generated wave loading. Two subjects were studied. First, the settlement data from the actual caisson breakwaters, measured over a time span greater than five years, were analyzed to determine the characteristics of the wave-induced settlement. The data obtained from the caissons were selected to build the west breakwater at the Jeju Outer Harbor. Secondly, a new numerical model consisting of two submodels: the wave and soil models, was proposed to simulate the wave-induced settlement of a caisson built on a sandy seabed. The dynamic wave load calculated by the wave model (2D-NIT) was used as the input data for the soil model (FLIP). This soil model can simulate both oscillatory and residual pore water pressures at the same time.

#### 2. Characteristics of the Soil and Waves at the Site

##### 2.1. Properties of the Soil and Waves at the Site

The results of a boring investigation conducted at the west breakwater construction site in the Jeju Outer Harbor indicated the seabed ground level was DL (datum level) (-) 12.78∼(-) 22.62 m. The bedding in this region consists of marine alluvial deposits, weathered bedrock deposits, and soft rock from the face of the seabed surface. The marine alluvial deposits, the main subject of this study, consist of fine to medium sand and silty sand, based on USCS (unified soil classification systems) guidelines. It has been transported by waves and distributed between 0.3 and 17.5 m, with an average of approximately 10 m. This sediment sand contains a significant amount of carbonate sand formed by fish bones and shellfish. The *N* value of the standard penetration test ranged from 9 to 25, which indicates a generally dense or very dense relative density. The variations were based on the location and the depth. Interestingly enough, according to the soil profile in Figure 2, no soft clay layer exists in this area.

A self-recording wave gauge was installed to observe the height and direction of the waves at the west breakwater construction site. The impact of the wave reflection, due to the structure of the caisson body, was minimized because the measuring point was 500 m from the front of the breakwater. This space between the gauge and the breakwater was large enough to eliminate the wave reflection effect. Table 1 summarizes the observation of waves over 5 m from 2002 to 2006. During the monitoring period, one mega typhoon, Maemi, occurred. Maemi produced waves with a maximum height of 11.38 m (max. significant wave height of 5.7 m) and a maximum wave period of 7.8 sec (max. significant wave period of 8.4 sec). Additionally, as shown in Table 1, after Typhoon Maemi, more than 14 times the average number of waves over 5 m in height passed this area. On Dec 21, 2005, a storm hit this area, producing waves with a height of 11.01 m, with a maximum wave period of 9.7 sec.

##### 2.2. Measurement of the Caisson Settlement Results and Discussion

A total of 84 caissons were installed in the west breakwater. Caisson #5 was the first to be installed, followed by caissons #1, 2, 3, and 4, which were installed to connect to caisson #5. They were installed in the eastbound direction in sequence, from caisson #6. The caissons installed in the west breakwater can be classified into four groups: nonperforated, curved slit, enlarged slit, and different shapes curved slit.

When constructing a caisson breakwater, the stability must be judged, which is why it is important to establish a plan to measure its settlement. In the case of the Jeju Outer Harbor, only level measurements were conducted. Jeju TBM (No. 7) was used as a reference point. The caisson’s settlement was recorded after it was installed at the reference point. The settlement was measured at four edges, and those edges were numbered from one to four—starting from the seaward edge of the caisson.

Of the 84 total caissons installed in the west breakwater, Figure 3 shows the selected caissons’ (#4, 5, 7, 8, and 9) settlement data. These selected caissons were installed from approximately 2003 through 2006 and have experienced wave actions. The vertical dotted lines in Figure 3, which are very important, represent the date of the occurrence of a high wave that could have affected the caisson’s settlement. As shown in Figure 3, the final settlement varied for each caisson, but the settlement trend was very similar across all the caissons. The initial immediate settlement occurred after the installation of the caisson due to its weight. Then, creep settlement occurred steadily due to the properties of the sand in this area. Then, on particular dates, such as September 12, 2003, a significant amount of settlement occurred rapidly. It is unlikely that this rapid, significant settlement was caused by the material properties of the sand because the settlement progressed greatly in a short period of time. Rather, this settlement may have been induced by wave actions. The settlement occurred at the same time that a large wave hit.

The settlement of caisson #8, which was installed just three months before the typhoon, increased from 26 cm to 42 cm after the typhoon. The settlement of caisson #9, which was installed just one month before the typhoon, increased dramatically from 14 cm to 45 cm after the typhoon. Other caissons also exhibited a significant amount of wave-induced settlement, though of a lesser magnitude than caisson #9, because they were installed much earlier.

#### 3. Present Numerical Model

##### 3.1. Wave Model (2D-NIT)

To determine the dynamic wave force (including breaking waves) on the seabed and the breakwater, the present study used a program called the 2D-NIT (two-dimensional numerical irregular wave tank) model [21,22]. The model used the solutions to the viscous and incompressible Navier-Stokes equations for a two-phase flow (water and air) model, and the volume of the fluid (VOF) method was used to treat the free surface of the water. The numerical wave tank that was adopted in the 2D-NIT model is shown in Figure 4. As for the open boundary conditions, fictitious dissipation zones (La) were added to the left and right sides of the computational domain. Their thickness is 2*L*, with *L* being the wavelength, in order to absorb the wave energy. The internal wavemaker is located in front of the fictitious dissipation zone on the left side of the computational domain.

The basic equation of the 2D-NIT model consists of the continuity equation (1) and Navier-Stokes’ equations of motions 2 and 3, which were expanded to a porous media by Sakakiyama and Kajima [23] described as follows:where *t* is time; *x* and *z* are the horizontal and vertical coordinates; *u* and are the horizontal and vertical components of flow velocity; *ρ* is the density of the fluid; *p* is the pressure; is the sum of the dynamic viscosity coefficient and the eddy viscosity coefficient; is the acceleration of gravity; is the volume porosity; *λ*_{x} and *λ*_{z} are the directional area porosity; *D*_{x} and *D*_{z} are the wave energy dissipation coefficients for directions *x* and *z*; *S*_{ρ}, *S*_{u}, and are the source terms to generate a wave in the computational domain.

The 2D-NIT program uses the VOF method to trace the free surface. In the VOF method, the interface between the water and the air phase is modeled according to the VOF function; i.e., the VOF method evolves the volume of water in each cell over time instead of directly tracking the free surface itself. The VOF function *F* is the volume fraction of a fluid that is a constant physical quantity and has a value of 0 ≤ *F* ≤ 1.

##### 3.2. Soil Model

The FLIP program, a finite element analysis program for the liquefaction process that is based on the effective stress analysis method Iai et al. [24], was adopted to analyze the storm wave-induced settlement behavior of the caisson breakwater. This program was originally developed for the dynamic soil response under earthquake loading. The effective stress model for these soils is the multiple shear mechanism model that was originally proposed by Towata and Ishihara [25]. As shown in Figure 5, this model is represented by a movable point located within the circular fixed boundary. It is defined in the shear strain space and connected to the boundary with an infinite number of virtual springs. Each spring corresponds to a virtual simple shear mechanism and has one of a number of various orientations. The relationship between the force and the displacement of each spring follows the hyperbolic-type load-displacement pattern. The displacement of the movable point from the center represents the mobilized shear strain, and the resulting forces acting on the point represent the shear stress that is induced in the soil.

The excess pore water pressure is generated as a function of the cumulative shear work [24]. Equation (3) considers the effect of positive dilatancy when taking the cyclic mobility behavior into account using the concept of a liquefaction front, as shown in Figure 6:in whichwhere *S*_{0} is a parameter to be defined by a function of shear work; *m*_{1} is the inclination of the failure line, defined by the shear resistance angle *ϕ*_{f}′ as *m*_{1} = sin*ϕ*_{f}′; *m*_{2} is the inclination of the phase transformation line, defined by the phase transformation angle *ϕ*_{p}′ as *m*_{1} = sin*ϕ*_{p}′; *m*_{3} = 0.67*m*_{2}. The auxiliary parameter, *m*_{3}, is introduced to ensure a smooth transition from one zone to the other and is determined as a balance of the smoothness and the realism of the stress path shape.

In Figure 6, *S* is a state variable (*S* = *σ*_{m}′*/σ*_{m0}′) under undrained conditions with a constant total confining pressure, and *r* is the shear stress ratio (*r* = *τ*/(−*σ*_{m0}′)). The initial effective mean stress and deviatoric stress are defined by *σ*_{m0}′ = (*σ*_{x0}′ *+* *σ*_{y0}′)/2 and *τ* = (*σ*_{1}′ − *σ*_{3}′)/2 = (*τ*^{2}_{xy} *+* ((*σ*_{x}′ − *σ*_{y}′)/2))^{0.5}, respectively. The model can simulate a rapid or gradual increase in the cyclic strain amplitude to the order of several percentage points under undrained cyclic loading. The program has been verified in many numerical simulations of structure damage induced by earthquakes and liquefaction [26–29].

##### 3.3. Verification of the Present Numerical Model

Verification of the present numerical model is a necessary process before applying it. To verify the present numerical model, two previous studies available in the literature, Mizutani et al. [9] and Jeng et al. [12] were used. Mizutani et al. [9] conducted a series of experiments for regular waves to investigate the wave-seabed-structure (submerged breakwater) interaction. The water-free surface elevation and the dynamic pore water pressure inside the sandy seabed foundation under the breakwater were recorded during experiments. Mizutani et al. [9] also developed a combined BEM-FEM model. Modified Navier-Stokes equations were used to solve the flow inside the porous media, and Biot’s equation was applied to solve the poroelastic media. Jeng et al. [12] developed the integrated model (PORO-WSSI II) for wave-seabed-structure (WSSI) by combining the Volume-Averaged Reynolds-Averaged Navier-Stokes (VARANS) equation for wave motions in a fluid domain and the porous media flows in structures, and the dynamic Biot’s equations for a porous elastic seabed. Jeng et al. [12] applied their model to investigate the dynamic response of a composite breakwater on a seabed.

The properties of the sand and breakwater provided by Mizutani et al. [9] and Jeng et al. [12] listed in Table 2 were used for simulation of the present numerical model. The soil model (FLIP) in the present model needs more parameters related to the liquefaction process. But in this verification, those parameters were assumed as zero because (i) sand has no fines, and (ii) Mizutani et al. [9] and Jeng et al. [12] only considered oscillatory pore water pressure in seabed based on Biot’s equation. The regular wave of height 3.0 cm and period of 1.4 sec in a still water depth 30 cm above the sand bed provided by Mizutani et al. [9] and Jeng et al. [12] were used for simulation.

Figure 7 shows a schematic sketch of the numerical wave-soil tank used for verification of the present model. Four locations for the value of *η*/*H* (*η* = water surface elevation above SWL and *H* = wave height) were selected and named as “a, b, c, and d.” Location “a” is in front of the breakwater, “d” is behind the breakwater, and “b” and “c” are edges of the breakwater crown. Four locations for the pore water pressure selected and named as *A*, *B*, *C*, and *D*. Location “*A*” was used for recording the pore water pressure in the middle of the submerged breakwater and “*B*, *C*, and *D*” for recording the pore water pressure inside the seabed in front of, middle of, and behind the breakwater.

Figure 8 shows the computed dimensionless water surface elevation (the value of *η*/*H*) with the 2D-NIT model in the present model together with Mizutani et al.’s experiment [9] and Jeng et al.’s numerical result [12]. The circle, blue line, and red line indicated Mizutani’s experiment, Jeng’s numerical result, and the computed one with the 2D-NIT model, respectively. As illustrated in Figure 8, the agreements at “a” and “b” are good, while little differences are observed at “c” and “d” that are located behind the breakwater in between the computed and Mizutani’s experiment values. This may be accounted to the boundary condition at the seabed. In the experiment of Mizutani, it is assumed as permeable, whereas in the numerical analysis, it is assumed as an impermeable rigid surface. Anyway, the difference is small enough, and the wave model used in this study can successfully simulate wave propagation in a viscous fluid.

**(a)**

**(b)**

**(c)**

**(d)**

Figure 9 shows the computed dimensionless pore water pressure (*p*_{s}/*ρgH*, in which *p*_{s} = pore water pressure and *ρ* = density) with the FLIP in the present model together with Mizutani et al.’s experiment in Mizutani et al. [9] and Jeng et al.’s numerical result [12]. The circle, blue line, and red line indicated Mizutani et al.’s experiment, Jeng et al.’s numerical result, and the computed one with the 2D-NIT model, respectively. The numerical perdition of the pore water pressure overall agrees well with the results of Mizutani’s and Jeng’s works. From this result, it can be stated that the present soil model can simulate the wave-induced pressure breakwater and its seabed accurately.

**(a)**

**(b)**

**(c)**

**(d)**

##### 3.4. Application of the Present Numerical Model to Caisson Breakwater

###### 3.4.1. Modeling of the Target Breakwater

The caisson breakwater shown in Figure 10 was modeled for the numerical analysis, as shown in Figure 11. The rubble mound is 6 m in height and 72 m in width with a slope of 1 : 1.5. The caisson is 26 m in height and 24 m in width. The thickness of the wave dissipation block (tripod) is about 1∼2 m. The sandy soil layer is with 16 m thick. It is placed as close to the target site ground (Figure 2) as possible, and it is located under the rubble mound. A layer of weathered rock is placed below the sandy layer.

**(a)**

**(b)**

Figure 12 shows the modeled meshes of the composite breakwater and the seabed. The seabed underneath the caisson breakwater is modeled with a smaller mesh. Figure 13 indicates the results of the mesh convergence test. It is the relationship between displacements of the caisson breakwater and the number of elements. The displacement was increased with the number of elements and they were converged when the number of elements was more than approximately 2,300. Therefore, 2,326 elements were used to make the mesh in this study.

Dynamic wave pressures acting at the boundary of the caisson breakwater and the seabed are calculated at all points. The dynamic responses of the caisson breakwater and the seabed are also estimated at all locations. As a matter of convenience, several points were selected to represent the general behavior of the caisson breakwater structure and the seabed. In Figure 12, *N*1 and *N*2, at the top of the caisson, mark the locations of caisson settlement. Element numbers 1∼30 indicate the representative wave pressure output points. Elements *E*1 ∼ *E*4 indicate the representative output elements of seabed behavior.

###### 3.4.2. Wave Pressure Estimation Using 2D-NIT

The northern direction of an irregular wave with a maximum height of 7 m and a maximum period of 11.0 sec was used as the source wave in this analysis. This source wave was selected after the characteristics of Typhoon Maemi (Table 1) were considered. A water depth of 13.3 m was used. This water depth was determined by integrating the tide level for the duration of the wave in this area and then adding the depth to the mean.

The numerical wave channel used in the present study is 1,750 m long and 60 m tall. Fictitious dissipation zones of 500 m are located to the left and right sides of the computational domain in order to absorb the wave energy. The lattice distance *Δx* is 0.5 m at the section where the breakwater was installed and 1.0 m at the other sections, while *Δz* chooses to be 0.5 m. The interval of time is Δ*t* = 0.01∼0.20 sec. The total calculation time is set at 600 sec. This given time is ample enough to study how the seabed’s behavior changes in response to the wave.

Dynamic wave pressures were obtained for all of the points on the breakwater structure and its surrounding seabed. This pressure data was used for the input data (the external force) in the seabed analysis program (FLIP). Figure 14 shows a typical wave load history: (a) at point No. 3 on the caisson side above the water level, (b) at point No. 10 on the caisson side below the water level, (c) at point No. 16 on the rubble mound, and (d) at point No. 25 on the seabed.

**(a)**

**(b)**

**(c)**

**(d)**

As seen in Figure 14, a significant number of nonuniform dynamic wave pressures are repeatedly acting on the breakwater except at the part over the hydrostatic water level, such as point No. 3. Additionally, a significant number of dynamic wave pressures are acting nonuniformly on the seabed (point No. 25) in front of the breakwater. These wave pressure distributions reveal some very important information. Some cases of breakwater design only consider the dynamic wave pressure on the breakwater structure itself and assume the wave pressure is uniform. There is a considerable difference between making assumptions about the wave pressure and actually studying it. Thus, the present study considered the real dynamic wave pressure acting on the seabed and the breakwater structure.

###### 3.4.3. Seabed Behavior Analysis Using FLIP

Before the dynamic wave loading response analysis was performed, static analysis was conducted where gravity was used to simulate the initial stresses acting on the seabed. The wave pressures obtained by 2D-NIT were assigned at the boundary of the domain as the input motion. The shear stress induced by the wave velocity on the surface of the seabed was also incorporated as part of the input source. The seabed was analyzed in an undrained condition. The peak period of an earthquake may vary over the range 0.1s∼1s smaller than the range of wave period. Owing to the smaller periods in the earthquakes, pore water pressure accumulates faster than in the case of waves. Mutlu Sumer [30] mentioned that the cyclic shear stresses cause the soil to undergo cyclic shear deformations, in exactly the same way as in the case of waves leading to pressure buildup, if the soil is undrained. It was assumed that the vertical and horizontal displacements were fixed on the bottom boundary. The vertical displacement was only allowed on the side boundary. The material parameters used for the analysis are shown in Table 3. They were determined based on field and laboratory test results.

*(1) Excess Pore Pressure Ratio*. The excess pore pressure ratio is the ratio between the initial effective stress and the excess pore pressure (hereinafter EPP). When the EPP ratio reaches up to 1, and we can define the soil element, it becomes liquefied. Figure 15 shows the EPP ratio in seabed elements *E*1, *E*2, *E*3, and *E*4. The EPP ratio of all elements increases with time. This is especially true for element *E*3, as it shows a relatively high increase in its EPP ratio. The EPP ratio of element *E*3 increases up to 0.9, while the other elements are around 0.5. These values are not 1, so the sand layer may not liquefy completely. However, it is close to reaching a value of 1, so we can expect the effective strength of the sand layer to decrease without completely vanishing. This is called “partial liquefaction” [18,19].

The buildup of EPP on the seabed and underneath the caisson may be caused by the combination of both the wave and the caisson’s movements [17,20]. Oumeraci’s [17] study suggested an especially interesting result for the pore pressure of the seabed where the breakwater was placed during the wave tank model that was undertaken in the condition of wave height *H* = 0.4 m, wave period *T* = 6.5 sec, and water depth *hs* = 1.6 m. The behavior of the seabed ground layer is the same in accordance with the behavior of the wave, crest, and trough. However, the wave does not directly trigger the excess pore pressure in the lower ground area beneath the caisson. Alternately, the behavior of the seabed underneath the caisson breakwater was essentially governed by the wave-induced caisson motion. If the caisson moves upward, the ground shows a negative oscillatory excess pore pressure, while it shows a positive oscillatory excess pore pressure in cases where the caisson moves downward. During the caisson’s movements, due to the wave, the excess pore pressures in the ground underneath the caisson keep increasing.

Figure 16 shows the final state (after 600 seconds of wave motion) of the EPP distribution. The critical part of the EPP ratio is the edge of the rubble mound on the seaward and shoreward sides. This may be related to the location, which is on the border of the seabed and the breakwater structure. In this area, the EPP had significant buildup due to the wave and the breakwater motions. Generally, near the edges, the shear stress can be much higher.

*(2) Effective Stress Path*. Figure 17 shows the effective stress path in the *p*′ − *q*′ diagram. The effective stress path is a very useful way to track the stress changes in the ground. During a wave, the downward propagation of the wave force through the seabed generates shear stresses and strains that are cyclic in nature [31]. If a cohesionless soil is saturated, excess pore pressures may accumulate during wave shearing, and the effective stress path moves toward failure. As shown in Figure 16, as the EPP increases, the effective stress path moves toward the Mohr-Coulomb failure line for all of the soil elements. Element E3 in the seabed, related to the EPP ratio shown in Figure 15, has a relatively high increase in its EPP ratio and is closer to the failure line than any other element. A higher increase in the EPP results in a higher decrease in effective stress.

*(3) Settlement of Caisson Breakwater*. Figure 18 shows the settlement deformation of the caisson breakwater and the area around the seabed at the end of the analysis, after 600 seconds. For a better view, this figure was magnified twenty times. The settlement deformation mainly occurs in the top soil layer beneath the caisson breakwater. This settlement was caused by the effective stress decrease in the soil layer, which stems from the EPP increase. When a sand layer is subjected to cycles of shear strains under undrained conditions, excess pore pressure may be generated in each load cycle leading to softening and accumulated deformation [32]. An increase in the pore pressures leads to a transfer of stress from the soil skeleton to the pore pressure. This precipitates a decrease in the effective stress and shear resistance of the soil. If the shear resistance of the soil becomes less than the driving shear stress, the soil can undergo large deformations [31].

In this numerical simulation, the EPP ratio increases and the effective stress decreases in the seabed sand layer when it is subjected to cyclic shearing stress induced by the wave force (Figures 15 and 17). A decrease in the wave-induced effective stress is a result of the seabed ground’s tendency to decrease in volume (settlement). Castro [33] classified settlement deformation into two categories depending on the shear stress and shear strength: (1) when the driving shear loads are greater than the residual shear strength of a liquefied soil deposit, a loss of stability can result in extensive ground failures or flow slides, and (2) when the driving shear stress is less than the residual shear strength, limited shear distortion and settlement occur without the soil mass losing its stability. This present study is corresponding to the latter because the EPP values do not reach 1, so the sand layer is not completely liquefied.

Figure 19 shows the combined actions of the wave pressure on element No. 10 in the caisson (above) and the settlement behavior of the caisson (below) over time. From this figure, evidence of the effect of a wave on the settlement of the caisson breakwater is clearly visible. The settlement of the caisson has not happened at the initial start time (*t* < 40 sec) because there was no wave pressure acting on the caisson. As time goes on, the wave pressure starts acting on the caisson wall and induces the settlement of the caisson due to the increase in the excess pore pressure on the seabed. The settlement continuously accumulates under the wave pressure and at around 340 sec, the time where the maximum amount of wave pressure is acting on the caisson; the maximum settlement has occurred. After that, the settlement is less significant because the pressure is less significant.

*N*1 and *N*2 are points on the top of the caisson. As seen in Figure 19, the trends of settlement at *N*1 and *N*2 are similar. The settlement of *N*2 (landward) is bigger than that of *N*1 (seaward). Figure 20 shows the settlement measured at four edges of caisson No. 6 that was installed in the west breakwater of the Jeju Outer Harbor. Coincidentally, the settlements at the edges of the landward side (numbered as 3 and 4) are bigger than those at the edges of the seaward side (numbered as 1 and 2). This may be related to the stress concentration at the landward seabed due to the high wave pressure acting on the seaward side of the caisson.

#### 4. Conclusion

This study was conducted to characterize wave-induced settlement in coastal structures (including the ground underneath these structures) by analyzing the settlement data measured in actual breakwaters over a long period of time (more than 5 years). To determine the effect that the waves have on seabed settlement, the wave history that was recorded during the construction period was also analyzed. These analyses show that waves clearly influence the settlement of caisson breakwaters. Wave-induced settlement in the breakwaters becomes very significant and rapid when large waves attack, such as during a typhoon.

A new numerical model was also suggested to simulate the storm wave-induced seabed settlement underneath the caisson breakwater, which could be simulated qualitatively. The new model combines the wave model (2D-NIT) and the soil model (FLIP). The dynamic wave load calculated by the 2D-NIT was used as the input data for the soil model. The model can simulate the oscillatory and residual pore pressures increase and the effective stress decrease in the seabed sand layer when the layer is subjected to cyclic shearing stress induced by wave force. The decrease in the wave-induced effective stress results from the seabed ground’s tendency to decrease in volume (settlement). The numerical analysis and field measured data both show clear evidence of the effect of waves on the settlement of the caisson breakwater.

#### Data Availability

The research data used to support the findings of this study were supplied by 2D-NIT and FLIP under license and so cannot be made freely available. Requests for access to these data should be made to the corresponding author.

#### Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this article.

#### Acknowledgments

This research was supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (No. 2020R1I1A3067248).