1. Introduction
Wave runup around offshore structures in extreme conditions is of practical significance due to the possible onset of wave slamming into the underside of a platform, or greenwater on deck phenomena for vessels. The latter, characterized by a large amount of water impacting onto the deck and/or deck structures, is an important issue for the efficient operation and, ultimately, the survivability of floating production, storage and offloading (FPSO) vessels (HSE 2005). This is a typical highly nonlinear wave–body interaction problem in which the nonlinearities could result from dynamic coupling among waves of different frequencies, waves and vessel motions and different modes of vessel motions. Accurate understanding and prediction of these fully nonlinear effects are still challenging, owing to the complexity of the hydrodynamics involved.
Nonlinear wave fields scattered from vertical surface-piercing cylinders have been the focus of much previous work using a perturbation approach, including Eatock Taylor & Huang (Reference Eatock Taylor and Huang1997) and Kim & Yue (Reference Kim and Yue1990) for nonlinear wave diffraction up to second order (Faltinsen, Newman & Vinje (Reference Faltinsen, Newman and Vinje1995) and Molin et al. (Reference Molin1995) captured third harmonic forces but not free surface elevations). In addition, fully nonlinear wave–cylinder interactions have been investigated numerically e.g. in Bai & Eatock Taylor (Reference Bai and Eatock Taylor2007), Chen et al. (Reference Chen, Zang, Hillis, Morgan and Plummer2014) and Jian et al. (Reference Jian, Cao, Lo, Huang, Chen, Cheng, Gu and Li2017), and experimentally e.g. in Riise et al. (Reference Riise, Grue, Jensen and Johannessen2018a), Riise et al. (Reference Riise, Grue, Jensen and Johannessen2018b) and Chen et al. (Reference Chen, Zang, Taylor, Sun, Morgan, Grice, Orszaghova and M Tello2018). Second- and third-order diffraction are found to be important in certain applications (i.e. for certain structure sizes and within certain flow regimes).
In contrast to vertical cylinders, the contribution from higher-order nonlinear components of wave–structure coupling for a ship-shaped geometry has remained unclear. Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006) showed that the measured runup on the bow of a simplified fixed model with head-on waves was entirely consistent with second-order diffraction theory. The sum and difference components were both important and the difference term corresponded to a substantial and persistent mound of water formed around the bow. Similar observations have been reported by Mai et al. (Reference Mai, Greaves, Raby and Taylor2016) in their experiments for fixed FPSO-shaped bodies.
For a freely floating structure, Stansberg & Karlsen (Reference Stansberg and Karlsen2001) found that the largest relative wave–vessel motions in their dataset are close to Rayleigh distributed because the pitch motion is overpredicted, which cancels out the underprediction of the free surface elevation. Stansberg & Berget (Reference Stansberg and Berget2009) and Schiller et al. (Reference Schiller, Pâkozdi, Stansberg, Yuba and e Silva2014) concluded that the nonlinearities in the relative wave kinematics mainly result from the nonlinear asymmetry of the incident wave elevation and the nonlinearity associated with local waves diffracted/radiated from a floating FPSO. The vertical vessel motion (mainly due to the pitch motion), however, is found to be well predicted by linear analysis, except for waves that are clearly longer than the vessel. Buchner (Reference Buchner2002) argued that the pitch vessel motion is nonlinear due to the effect of water on deck, and has an effect opposite to the nonlinearity in local wave peaks that leads to an increase in the relative wave–vessel motions. In contrast, Ruggeri et al. (Reference Ruggeri, Watai, de Mello, Sampaio, Simos and e Silva2015) stated that the vessel motion can be nonlinear regardless of the water on deck, in particular in steep waves. It is worth noting that greenwater overtopping was considered in these references. Obviously, a more thorough investigation is required to clarify the influence of nonlinear effects.
Our present work aims at improving understanding of the underlying physical features of the nonlinear solutions, making the development of a simpler model for predictions of extreme wave runup viable. Such a method should be sufficiently accurate in terms of the input to hydrodynamic load calculations for complex greenwater events. To achieve this, a particular focus is placed on examining if second-order theory is accurate enough for the prediction of wave diffraction and radiation by a freely floating ship-shaped body, with emphasis on accurate extraction of the harmonics and the contributions from various underlying physical processes.
It is noted that such an engineering-type method is attempted by Stansberg et al. (Reference Stansberg, Berget, Hellan, Hermundstad, Hoff, Kristiansen and Hansen2004) and Schiller et al. (Reference Schiller, Pâkozdi, Stansberg, Yuba and e Silva2014) in which an empirical 
$10\,\%$ increase (determined from their model test results) is added to linear diffraction/radiation solutions. The nonlinearity in the incoming wave is, however, simulated by a second-order random wave model. A reasonably good agreement to the experimental measurements is obtained for the statistics. Unlike Stansberg et al. (Reference Stansberg, Berget, Hellan, Hermundstad, Hoff, Kristiansen and Hansen2004) and Schiller et al. (Reference Schiller, Pâkozdi, Stansberg, Yuba and e Silva2014), among others, a wave-by-wave analysis, rather than a statistical analysis on nonlinear random waves, is carried out in this work by using transient wave groups as input waves to excite nonlinear wave scattering from a simplified FPSO hull with a vertical bow. The FPSO is either fixed or freely floating.
The broadbanded structure of ocean waves is retained using a wave group. This also makes computationally expensive numerical simulations tractable and minimizes the contamination by reflected waves from the sides of finite wave tanks when compared to results for regular waves and prolonged random sea experiments. Most attractively, the use of transient wave groups works very well in tandem with the phase separation method (Fitzgerald et al. Reference Fitzgerald, Taylor, Eatock Taylor, Grice and Zang2014), allowing for the extraction of individual nonlinear harmonics. The nonlinear contributions from each nonlinear process described above can then be examined in detail.
The phase separation method assumes that the nonlinear process of interest can be described by a Stokes-type perturbation expansion in both frequency and wave steepness. That is, the magnitude of each harmonic at frequencies that are close to integer multiples of the fundamental frequency is proportional to the linear component amplitude (regarded as a slowly time varying envelope of linear components when generalized to wave groups) to the corresponding power. The form associated with each harmonic is obtained from a linear combination of the time histories resulting from inputs with constant increments in phase, and their Hilbert transforms (Fitzgerald et al. Reference Fitzgerald, Taylor, Eatock Taylor, Grice and Zang2014). In the case where four wave groups are used, these have the same linear envelope but with each Fourier component shifted by a relative phase of 
$0^{\circ }$, 
$90^{\circ }$, 
$180^{\circ }$ or 
$270^{\circ }$, respectively. It has already been demonstrated that such generalized Stokes expansions are appropriate for analysing nonlinear loading on a bottom mounted cylinder in water waves of moderate steepness (Fitzgerald et al. Reference Fitzgerald, Taylor, Eatock Taylor, Grice and Zang2014) and in locally violent nearly breaking waves (Chen et al. Reference Chen, Zang, Taylor, Sun, Morgan, Grice, Orszaghova and M Tello2018), as well as for the resonant response in the gap between two parallel fixed boxes (Zhao et al. Reference Zhao, Wolgamot, Taylor and EATOCK TAYLOR2017). The generalization of this Stokes-like approximation for extreme wave runup will be tested against well-resolved computational fluid dynamics (CFD) simulations and diffraction analysis up to second order in this study. The application of the methodology to the FPSO problem is further complicated by the vessel motion and its interaction with nonlinear waves.
2. Numerical and experimental analyses
Both physical and numerical experiments are performed to investigate the extreme wave runup around an FPSO model, although this work puts emphasis on the numerical analysis. In addition to the measurements at certain key locations, we also examine the global free surface motion across the full spatial domain. The latter can be achieved more easily by running numerical simulations. However, the physical interpretation is equally appropriate for the wave tank tests. Verification of the numerical results is carried out by comparisons to the experimental measurements. In this section, numerical and experimental details of primary importance are presented.
The experiments were carried out in a 4 m wide wave flume at Dalian University of Technology, China, with a constant water depth of 1 m. A simplified FPSO form with a semi-circular bow and stern (in plan) was used and was placed at a location approximately 17 m away from the wave paddle. The wave flume has a length of 60 m. In order to match the box-like shape of many real FPSOs, the length of the straight portion of the hull was set to 1.31 m, so the entire model was 1.6 m at 
$1:200$th laboratory scale. The floating model was adjusted to achieve the desired draft of 0.1 m. The width and the depth of the model FPSO were 0.29 m and 0.175 m, respectively. The mass of the floating model was evenly distributed (i.e. radii of gyration in roll and pitch were 0.0936 m and 0.4473 m, respectively). A soft mooring system was used, which consisted of 4 taut wires orientated diagonally outwards from the model to points on the tank sides; each wire was attached to a spring, and the spring to the tank edge. This system provided virtually no restraint for vertical motion (so heave, pitch and roll) but suppressed large motions in the horizontal plane (so surge, sway and yaw). The experimental set-up is shown in figure 1. A stiff plastic sheet was wrapped around the model hull to prevent wave overtopping onto the deck. Relative wave runup was measured at the various wave gauges located around (and fixed to) the hull; the location of the wave gauges are shown in figure 1(b).
Figure 1. Experimental set-up. (a) Photo from the experiment; (b) schematic overview of the experiment. WG stands for wave gauge – these were fixed to the hull, thus measuring relative motion. Additional free-field wave gauges (black dots) were fixed to the flume.
Only long-crested waves can be created in the flume, and the model was oriented parallel and at an angle of 
$30^{\circ }$ to the flume walls (i.e. the vessel headings considered were 0 
$^{\circ }$ and 
$30^{\circ }$). The assumed underlying sea state had a Joint North Sea Wave Project (JONSWAP) spectral shape with a peak wave period of 1.05 s. This is in line with the least desirable sea states for extreme wave runup at the bow identified by a response-based analysis for the model geometry using linear diffraction theory (Chen et al. Reference Chen, Taylor, Draper, Wolgamot, Milne and Whelan2019c). The maximum wave runup is found to be a result of the vertical vessel motions (mainly pitch motion here) occurring out of phase with the disturbed local free surface elevations. This indicates that the maximum relative wave runup does not necessarily coincide with the largest wave. Coupling of the incident waves with linear transfer functions (i.e. response amplitude operators (RAOs) of the relative wave runup) is of primary importance as well as the environmental condition itself. It is noted that waves with a period of 1.05 s would excite resonance in the pitch response of the vessel considered.
A moderate peak wave amplitude of 0.05 m was selected so as to avoid a possible extreme situation in which the vessel bottom would be lifted above the water surface, with possible bottom slamming on re-entry. Recall that the draft is 0.1 m in this work.
The incident wave group was calibrated to have a certain shape in time at the location of the model so as to excite the desired response in a given sea state (here with the spectral peak wave period 
$T_p = 1.05\ \textrm {s}$ and the significant wave height 
$H_s = 0.05\ \textrm {m}$ in laboratory scale, corresponding to 14.84 s and 10 m in full scale). In this work, the maximum linear relative vertical wave runup at the bow is desired (Chen et al. Reference Chen, Taylor, Draper, Wolgamot, Milne and Whelan2019c). The calibration is achieved in the experiments by applying the phase separation method in tandem with the iterative technique used in Chen et al. (Reference Chen, Stagonas, Santo, Buldakov, Simons, Taylor and Zang2019a). Utilizing the four-phase realizations (figure 2a), the linear free surface elevation (without the model in place) can be separated via addition and subtraction of the four signals and their Hilbert transform (Fitzgerald et al. Reference Fitzgerald, Taylor, Eatock Taylor, Grice and Zang2014). The paddle signal is iteratively corrected by comparing this linear component to the target/desired wave that has the desired shape as described above. The final linearized undisturbed wave shape is shown in figure 2(b) (black dotted line).
Figure 2. Time histories of the free surface elevation at the location of the model without the model in place, i.e. undisturbed waves. (a) The experimental time histories for all four full signals that have the same linear envelope but different global phases. (b) The corresponding linear free surface elevation extracted by applying the phase separation method to the full signals shown in (a). The dashed line indicates the time instant when the maximum runup occurred at the bow (at the left in figure 1).
The subsequent relative free surface elevations measured at the bow (i.e. wave gauge WG1 in figure 1) and mid-ship (i.e. wave gauge WG4 in figure 1) of the vessel are shown in figures 3 and 4, respectively (black dotted lines). Both full (top) and corresponding linear (bottom) results are shown. The maximum linear relative wave runup at the bow is 
$\sim$ 0.085 m. This runup would be excited by an equivalent severe sea state of 
$H_s \sim 0.085\ \textrm {m}$ and 
$T_p \sim 1.05\ \textrm {s}$ at laboratory scale were the vessel motion not considered (i.e. the response is measured relative to the position of the undisturbed water surface).
Figure 3. The relative vertical wave runup at the bow when the model was installed. (a) The total relative vertical wave runup and (b) the corresponding linear component. The vessel heading is 
$0^{\circ }$. Only pitch and heave motions are released in the CFD simulations though all 6 degrees of freedom (DoF) motions are allowed in experiments. The dashed line indicates the time instant when the maximum runup occurred at the bow (at the left in figure 1).
Figure 4. Same plot as figure 3 but for a location amidship of the vessel, WG4.
Using a 
$1:200$ length scaling from the tank to full field scale, this hypothetical underlying sea state has 
$H_s \sim 17\ \textrm {m}$ and 
$T_p \sim 14.84\ \textrm {s}$ at field scale. This corresponds to a 1 in 100 year return period sea state (Chen et al. Reference Chen, Taylor, Draper, Wolgamot, Milne and Whelan2019c). So this is a very severe and very steep sea state. If the phase separation method works for our analysis, it should be valid for many sea states of practical interest, at least for the simplified hull geometry used here. The possible application of the methodology to a wider range of wave conditions and vessel headings is being assessed, and will be reported elsewhere.
The experiment is reproduced using CFD-based simulations. Here, the OpenFOAM CFD scheme is utilized with the toolbox ‘waves2Foam’ for generating the incident wave groups (Jacobsen, Fuhrman & Fredsøe Reference Jacobsen, Fuhrman and Fredsøe2012). The governing Navier–Stokes equations are solved using a finite volume method. The water motion is captured by the volume-of-fluid (VOF) method. The experimental incident wave group is re-created with reasonable accuracy in the numerical wave tank (NWT) using an iterative method (see figure 2 and Chen et al. (Reference Chen, Stagonas, Santo, Buldakov, Simons, Taylor and Zang2019a) for details). We note here that the total length of the NWT is minimized to save computational cost, while the width and the depth are kept the same as in the experiments. The total length is set as 
${\sim }10\lambda _p$ with the model being placed 
${\sim }3.5\lambda _p$ away from the input boundary to allow wave propagation (Chen et al. Reference Chen, Stagonas, Santo, Buldakov, Simons, Taylor and Zang2019a). Here, 
$\lambda _p$ is the peak wavelength (1.72 m in this work).
The laminar flow model of OpenFOAM is used here as the wave runup at a vertical bow is assumed to be dominated by potential flow effects. This assumption is consistent with the excellent match between tank tests and potential flow modelling (second-order diffraction theory) obtained by Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006) for waves at the bow of a comparable fixed model. The maximum Reynolds number 
$Re$ 
$(=\omega {\eta ^2_m}/\nu$) here is approximately 
$5.98\times 10^4$, where 
$\omega$ is the peak wave angular frequency, and 
$\nu$ the kinematic viscosity; 
${\eta }_m$ is the maximum free surface elevation relative to the bow which is 
${\sim }0.1\ \textrm {m}$ in this study. Jensen, Sumer & Fredsøe (Reference Jensen, Sumer and Fredsøe1989) suggested that the laminar-to-turbulent transition would first occur at 
$Re \sim 10^5$ for oscillatory flow over a smooth wall. Hence, the use of laminar flow model is plausible here, especially when the oscillatory flow is transient. The generally good agreement between the numerical and experimental results as well as the diffraction solutions shown in the following sections justifies our assumption that the influences of viscosity, flow turbulence and skin friction on the main physical variables of interest are small. Likewise, the effect of surface tension is not numerically represented, since the Bond number (here 
$Bo =gL^2\rho /T_s = 3.4\times 10^5$ where 
$L$ is the length of the FPSO, and 
$T_s = 0.073\ \textrm {N}\,\textrm {m}^{-1}$ is the surface tension) for this problem is well above the range where surface tension is important according to Faltinsen & Timokha (Reference Faltinsen and Timokha2009). A static contact angle of 
$90^{\circ }$ is thus used to avoid modelling the surface tension forces between the wall and the fluid.
The whole computational domain is divided into several areas (in top view) and layers (in side view) with different levels of grid resolution. The actual mesh size used is determined by convergence tests; the mesh size along the wave propagation direction is chosen to achieve 
${\sim }170$ cells per peak wavelength to resolve propagating incident waves, and the vertical cell height is chosen to be 20 cells per peak wave height near the free surface to capture wave motions. The horizontal cell size is reduced to 
${\sim }340$ cells per wavelength in regions around the structure, and the cells are graded in the vertical direction so that the cells at the bottom and top of the tank are four times larger than those around the free surface. The dimensionless wall distance (
$y+$) is thus smaller than 5 based on the flat-plate boundary layer theory. It is found that further refining the mesh inside the boundary layer has negligible influence on the wave runup at the vertical bow (see also Chen et al. Reference Chen, Stagonas, Santo, Buldakov, Simons, Taylor and Zang2019a). The time step of the simulations is chosen to be runtime adaptive with a maximum Courant number not exceeding 0.25 to ensure stability. The total cell number is 
${\sim }40$ million, and the simulations within the time frame of interest (
${\sim }10\ \textrm {s}$ here) require 
${\sim }120$ thousand CPU hours (computed with 720 cores on a Cray XC40 supercomputer with 2.6 GHz Intel Xeon E5-2690 v3 ‘Haswell’ CPUs) for simulations with the floating model. For the same structure fixed, the computational cost with the same mesh resolution is reduced significantly to the order of several thousand core hours.
The calculated full and linear relative vertical wave runups at the bow and mid-ship of the vessel are compared with the experimental results in figures 3 and 4. The comparisons in terms of vessel motions are shown in figure 5. Satisfactory agreement is achieved, indicating that the numerical tool can provide accurate extreme runup around a FPSO as well as the vessel motions (The comparison has also been made for the fixed structure, with comparable agreement, but is not shown here for brevity). Therefore, it should permit a detailed exploration of the nonlinear scattering around the FPSO model. We note here that only pitch and heave motions are released in the simulations, although all 6 DoF motions were allowed in the experiments. This is reasonable as, for head seas, only pitch and heave body motions are important if the vertical runup at the bow is of interest.
Figure 5. The linear pitch (a) and heave (b) motions of the vessel. The total vessel motions are not shown, since we found that linear vessel motions actually dominate, as shown later in § 3. Negative pitch motion at 
$t = 66\ \textrm {s}$ means that the vessel bow is pitched down at this instant.
3. Nonlinear scattering around a FPSO
This section investigates the nonlinear wave fields diffracted and radiated from a fixed and a freely floating simplified FPSO model using the validated numerical scheme. It is noted that the validation for cases with a fixed body (by comparing to experiments) is not shown here, since this has been widely carried out elsewhere, such as Chen et al. (Reference Chen, Zang, Hillis, Morgan and Plummer2014) and Chen et al. (Reference Chen, Taylor, Draper and Wolgamot2019b). Nevertheless, the capability of the numerical scheme in simulating wave–fixed structure interaction is further confirmed by comparing to the diffraction solutions up to second order, as presented later in this section.
The phase separation method mentioned above will be adopted to extract the harmonic structure of the nonlinear wave fields across the full domain. Linear and second-order diffraction analyses are also carried out to verify the decomposition process, testing if the Stokes-like approximation is still appropriate even for this highly nonlinear system. It is noted that there are cross-terms in each harmonic component extracted using the phase separation method. The cross-terms have the same frequency, but a different (higher-order) dependence on the wave amplitude. For example, a third-order interaction of three linear components results in a term that scales as the cube of the linear wave amplitude, but has a frequency centred on the linear range. In general, all such subharmonic terms are assumed to be negligible. The generally good agreement between the diffraction analysis up to second order and the extracted CFD results shown in the following sections justifies this assumption for nonlinear waves scattered from a ship-shaped FPSO. We note that low-frequency subharmonic third-order terms may be of significance (i.e. can be of comparable magnitude to the second-order subharmonic response) for a catenary moored floater system, as presented in Orszaghova et al. (Reference Orszaghova, Taylor, Wolgamot, Madsen, Pegalajar-Jurado and Bredmose2020). Such systems typically have larger vertical natural periods and relatively weak hydrodynamic damping, i.e. the systems are highly resonant. But this is not the case in this work where the soft mooring system provided virtually no constraints on the vertical vessel motions.
The superscripts used in this work define the orders of harmonics and the potential flow solution, with (1) corresponding to the linear component, (2) the second harmonic etc.
3.1. Nonlinear diffraction from a fixed FPSO model
As mentioned in the Introduction, the four-phase runs allowed higher harmonics to be extracted using the phase separation method of Fitzgerald et al. (Reference Fitzgerald, Taylor, Eatock Taylor, Grice and Zang2014), with the results presented in figure 6 for the wave group incident from 
$0^{\circ }$. The disturbed free surface time history is analysed at the forward point of the bow. Both absolute dimensional and normalized results are shown in the time histories. The extracted harmonics are non-dimensionalized (right vertical axis) by dividing by the harmonic amplitudes associated with the Stokes expansion of incident bound waves on deep water, i.e. 
$\eta _0^{(1)} = A$, 
$\eta _0^{(2)} = A(kA)/2$, 
$\eta _0^{(3)} = 3(A(kA)^{2})/8$ and 
$\eta _0^{(4)} = A(kA)^{3}/3$. We take 
$A$ as the peak of the linear envelope of the incident wave group and 
$k$ is calculated based on the spectral peak wave period 
$T_p$.
Figure 6. Harmonic structures of the free surface elevation at the bow of the fixed FPSO in time (left), and their spatial distributions (right) at two typical instants (a,b) in order from top: total free surface elevation, low-frequency second harmonic component, linear (first harmonic) and high-frequency second harmonic component, third and fourth harmonics. The magnitudes of the second and higher harmonics are scaled up in the spatial plots to improve visibility; scale factors are 
$\times$15, 10, 50 and 125 for low-frequency second, high-frequency second, third and fourth harmonic components, respectively. The vessel heading is 
$0^{\circ }$. Red lines are CFD results and the dashed blue lines are diffraction solutions.
The diffraction solutions up to second order are also included using the numerical scheme outlined in Teng & Eatock Taylor (Reference Teng and Eatock Taylor1995). It can be seen that the agreement between the diffraction solutions and the corresponding harmonic components extracted from CFD simulations is in general good, with the values and phases matching very well. Although the comparison is highly acceptable, there are still some differences, especially for the low-frequency second harmonic component. This is not surprising as the wave generation in CFD (which is achieved by prescribing the velocity of the input boundary based on linear wave theory) will introduce some error waves that are not well separated from the low-frequency second harmonics due to the short length of the CFD domain. The cross-terms from a fourth-order interaction of four linear components, which scale as the fourth power of the linear wave amplitude but have a frequency range coinciding with low-frequency second harmonics, may also contaminate the results. However, this contribution is considered to be small, with the magnitudes being 
${\sim }100\times$ smaller. It should also be noted that the diffraction solutions are for a vessel in an infinite domain, while the CFD solution includes channel walls. The good agreement between the two methods for the linear and high-frequency second harmonic components suggests that the channel walls are unimportant for transient interactions with the bow. Overall, the Stokes-like expansion appears to be appropriate for this particular problem, as expected (so comparable to Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006)), allowing the exploration of the nonlinear diffracted wave field using the adopted separation method based on the CFD simulations.
The full interaction of the head-on wave group with the fixed FPSO is shown in supplementary movies 1–5 available at https://doi.org/10.1017/jfm.2020.1072. Snapshots of the interaction at two time instants (when the vertical runup is maximum at the bow and amidships on the side of the model, respectively) are shown in the field plots in figure 6.
There is relatively weak diffraction of the incoming linear component (the non-dimensional peak value is 
${\sim }1$). A strong and apparently localized interaction around the bow is clearly observed at the water surface from low-frequency second harmonic components, similar to that reported by Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006). This effect is driven by the products of linear terms (evident since the effect appears in the ‘quadratic’ terms in the diffraction solution), and decays rapidly in space away from the centre of the bow but lasts over a time scale comparable to that of the group rather than individual waves in the group. Note that the low-frequency second harmonic components under a unidirectional wave group lead to a set-down away from the model; here we have a local set-up at the bow. Both low- and high-frequency second harmonic components are important, with the contribution of second harmonics up to 
${\sim }37\,\%$ of the linear component (ratio of maximum values) at the bow, and up to 
${\sim }10\,\%$ amidships. All harmonics higher than second are small, less than 
${\sim }5\,\%$. Although these local scattered components are quite small in absolute terms, each is several times larger than the bound harmonics in the undisturbed incident waves. Similar phenomena have been observed by Mai et al. (Reference Mai, Greaves, Raby and Taylor2016) in their experiments, though their amplification in higher-order diffracted waves is relatively smaller. Their measurements were collected at a location slightly off the bow, while the location of the wave gauge in our work is at the stagnation point of the bow. It is worth noting that the waves were relatively longer in Mai's work (with wavelength/model length 
${\sim }2.3$ and 
${\sim }4.6$ compared to 
${\sim }1$ in our work).
Obviously, the phasing among the wave components is of importance to determine if higher harmonics actually increase or decrease the instantaneous total peak crest elevation. In this case, first harmonic scattering increases the linear crest elevation by 
$\sim$20 % (by comparing the peak values in linear local waves presented in the third plot (from top to bottom) of figure 6 left and linear undisturbed free surface elevation in figure 2b), and the second harmonic increases the crest value by a further 
$\sim$20 % (i.e. the contribution from both low- and high-frequency second harmonics), although they are not at their maximum values simultaneously. This is different to the results of Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006), where the linear and second-order terms are in phase for their crest focused wave groups, but the vessel model in Zang et al. (Reference Zang, Gibson, Taylor, Eatock Taylor and Swan2006) was larger relative to the incident wavelength. The incident wave group here is designed to maximize relative linear wave runup at the bow with the model free to move. It is not necessarily the optimum to maximize linear wave runup on the fixed model in the same sea state. This observation indicates that diffraction analysis up to second order is required to identify the form of the incoming wave within the group that would result in the nonlinear diffracted waves which maximize instantaneous vertical wave runup.
There is a striking similarity between the interactions with a vertical cylinder and with the simplified fixed FPSO here. A concentric wave field is formed at the bow associated with the run-up-wash-down process on the front face with the arrival of a wave crest. This was given the name Type-1 wave field in Sheikh & Swan (Reference Sheikh and Swan2005) and Swan & Sheikh (Reference Swan and Sheikh2015), where the interaction between steep waves and a surface-piercing column is detailed. The Type-1 wave radiates outwards mainly in the upstream direction. The second type of wave field, Type-2 wave, sees the increased water surface around the shoulder of the vessel where the circular front is joined to the rectangular hull. The induced disturbances propagate along the hull sides to produce a pair of symmetric, but non-concentric wave fronts (figure 6b), merging to form a notable mound of water at and behind the centre of the stern (cf. supplementary movies). The Type-2 wave is found to be more significant for higher harmonic shorter waves (figure 6b). If the incident waves are steep, a significant increase in the maximum crest elevation at the stern will be formed. This, in turn, may lead to greenwater at the stern where the freeboard of real FPSO vessels is usually smaller than at the bow. We note that the merging of the pairs of Type-2 waves only occurs downstream of the vessel here due to the relatively long hull. Merging upstream of the column (as presented in Swan & Sheikh (Reference Swan and Sheikh2015)) may break the harmonic analysis for a compact cylinder due to the violent interaction between the next incident wave crest and the Type-2 wave propagating in the opposite direction. But this upstream significant interaction does not seem to occur within the time frame of interest for the elongated hull shape and the transient wave group considered here, as indicated from the smooth full signal at the bow in figure 6 top plot.
Similar analysis has been carried out for the case with the incoming wave group incident from 30 
$^{\circ }$ relative to the vessel centreline; the results are shown in figure 7 and supplementary movies 6–10 (https://doi.org/10.1017/jfm.2020.1072). Similar conclusions to those in figure 6 can be drawn, although the scattered wave patterns are rather different. The approximation of the Stokes-type expansion works well, and the second harmonic component is significant while the harmonics higher than second are again relatively weak, and are perhaps negligible in the context of engineering models.
Figure 7. Same plots as figure 6 but for a fixed FPSO with incoming wave incident from 
$30^{\circ }$; (
$a_{S}$) and (
$b_{S}$) are the results for the diffracted wave fields only at the same time instant as (a,b), i.e. the undisturbed incident wave fields are subtracted from the total wave fields shown in (
$a_{T}$) and (
$b_{T}$), respectively. The magnitudes of the second and higher harmonics for both total and diffracted wave fields are scaled up in the spatial plots to improve visibility; for the colour code refer to figure 6.
Clear in the snapshots of the total wave fields at each harmonic in figures 7(
$a_{T}$) and 7(
$b_{T}$) are rings of higher harmonic waves scattered around the bow, and a significant interaction occurs off the weather side of the model as the higher harmonic bound components of the incident wave and part of the diffracted ones cross. There is strong diffraction around the stern of the model, leading to higher harmonic components propagating back upstream along the sheltered side (this is clearer for the relatively weak third and fourth harmonics; cf. supplementary movies). Again, the local scattered components at each nonlinear harmonic are several times larger than the bound components of the incident wave group. Corresponding results for the scattered and local generated wave fields only, so with the incident bound harmonics of the incoming wave group removed, are shown in figures 7(
$a_{S}$) and 7(
$b_{S}$), respectively. It can be seen that the structure of the linear diffracted wave fields is significantly different to that of the total linear wave fields. The opposite is true for higher harmonic diffracted wave fields, which are similar to the equivalent total fields, especially for third and fourth harmonics. This confirms that the linear incident waves dominate the linear components, while significant diffraction and harmonic wave generation occur for higher-order harmonics.
Obviously, the vessel heading has a significant effect on the local wave field, especially on the local wave fields around the sides and stern of the vessel. Detailed comparisons between the time series in figures 6 and 7 indicate that, for the forward point of the bow, the linear and fourth harmonic components are of comparable sizes for the two headings, while the behaviours of the high-frequency second and third harmonic components are rather different. The high-frequency second harmonic component decreases while the third harmonic increases as the heading increases from 
$0^{\circ }$ to 
$30^{\circ }$. This is different to the observation by Mai et al. (Reference Mai, Greaves, Raby and Taylor2016) that all harmonic components increase as the incident wave angle increases from 
$0^{\circ }$ to 
$20^{\circ }$ and the increase in the cases of 
$10^{\circ }$ is found to be the largest. Again, we note that the measurements were collected at a location slightly off the bow in their experiments while ours are at the forward points. Significant gradients in the interactions occur off the bow area as shown in the field plots of figures 6 and 7, hence, the observations of different trends for locations even a few centimetres away at laboratory scale are not surprising. This indicates that the exploration of the global behaviour across the full spatial domain is of importance.
3.2. Nonlinear diffraction and radiation from a floating FPSO model
The vessel motion is now released for the heave and pitch DoFs as these are most important for vertical runup, hence, greenwater overtopping at the bow. As for the local free surface elevation, the harmonic analysis is also carried out for the vessel motions; an example result for a wave group incident from 
$0^{\circ }$ is shown in figure 8. Both CFD and experimental results are included. It is found that the linear vessel (both heave and pitch) motion dominates, with contributions of higher harmonics to the ship motion smaller than 
${\sim }2\,\%$. This is also true for cases with the incoming wave incident from 
$30^{\circ }$ relative to the vessel centreline (harmonic structures not shown here for brevity). Therefore, comparisons in terms of linear components with the diffraction solutions are shown in figures 9 and 10 (experimental results from figure 8 are also included for reference), supporting the generalization of the methodology to freely floating soft-moored vessels. Note that pitch-down motion at the bow is defined as negative in this work.
Figure 8. Vessel motions for a floating FPSO with incoming wave incident from 
$0^{\circ }$. Here, 
$\theta$ and 
$z$ are pitch and heave motions, respectively. The superscripts indicate the order of harmonic components, and top plots are the total motions. Red and black dotted lines are extracted CFD and experimental results, respectively.
Figure 9. Linear vessel motions for a floating FPSO with incoming wave incident from 
$0^{\circ }$. Solutions from all three methods are included. Dashed lines indicate the time instant when the maximum runup occurred.
Figure 10. Same plots as figure 9 but for a floating FPSO with incoming wave incident from 
$30^{\circ }$.
Harmonic structures of the free surface elevation at the bow in time, and their spatial distributions at two typical instants, are shown in figures 11 and 12 for vessel headings of 
$0^{\circ }$ and 
$30^{\circ }$, respectively. The full interactions can be found in the corresponding supplementary movies 11–15 and 16–20 (https://cutt.ly/UhsHvSK), respectively. It can be seen that the maximum runup, i.e. maximum relative wave–vessel motion, occurs due to the vertical motion occurring out of phase with the local wave crests, as expected (dashed lines in figures 2–5, 9 and 10). It is worth mentioning that the time instant when the maximum runup occurs is shifted by 
$0.8{T_p}$ from 
$t = 65.16\ \textrm {s}$ for a fixed vessel to 
$t = 66\ \textrm {s}$ for a freely floating vessel. The former corresponds to the arrival of the maximum wave crest and the latter the passage of a second slightly smaller wave crest, as shown in figure 2. The maximum downwards vessel motion is excited by the former and is building up to its maximum at the latter time instant.
The diffraction solutions up to second order as well as the total local wave fields are also included (dashed blue lines). It can be seen that satisfactory agreement is obtained (as for fixed structures), which further confirms the validity of the adopted methodology in this work, i.e. the phase manipulation in tandem with CFD works well for extracting nonlinear higher harmonics.
It can be seen from field plots that the wave scattering patterns are similar. However, the amplification in local waves observed for the fixed vessel (figures 6 and 7) is weakened by the radiated waves associated with the vessel motions. A reduction of 
${\sim}20\,\%$ in both linear and second harmonics is found for the free surface elevation at the bow. In contrast, an increase is observed for the third and fourth harmonics, although their magnitudes are still relatively small. This indicates that the linear vessel motion actually dominates the vertical runups for freely floating vessels and the shorter radiated waves play a rather small role.
The harmonic structures of the relative vertical wave runup at the bow are shown in figure 13. Both experimental (measured by the wave gauge WG1 attached to the vessel, see figure 1) and CFD results (the combined results of the local waves and the vertical vessel motions, as shown in figure 14) are included. The diffraction solutions are not included for brevity, for which refer to figures 9–12. As with fixed structures and local waves as well as vessel motions, generally good agreement is obtained, as expected. A set of twin-wire resistive wave gauges, consisting of two parallel steel rods fixed at a set distance apart, were used in the experiments. The exact location (between the rods) at which the measurements are recorded is hence hard to determine. As mentioned above, the behaviour of the nonlinear wave fields (especially higher harmonic components) even centimetres away can be very different. Thus, slight differences in higher harmonics are expected. The differences observed for the case with a heading of 
$30^{\circ }$ may also result from the fact that the roll motion is released in the experiments, which is not represented in the CFD simulations. The roll motion may play a significant role for wave runup at a location off the centreline, however, the contribution of roll to the behaviour at the centre of the bow WG1 is found to be relatively small.
Figure 13. The harmonic structures of the relative vertical wave runup at the bow with incoming wave incident from 
$0^{\circ }$ (a) and 
$30^{\circ }$ (b). Red and black dotted lines are extracted CFD and experimental results, respectively. Dashed lines indicate the time instant when the maximum runup occurred.
Figure 14. The disturbed local wave (black lines; incident+scattered+wave making) and the vertical vessel motion, VM, due to heave (red lines) and pitch motions (blue lines). The results for the heading of 
$0^{\circ }$ and 
$30^{\circ }$ are shown on the top and bottom, respectively. Dashed lines indicate the time instant when the maximum runup occurred.
4. Conclusions
In this work both numerical simulations solving the Navier–Stokes equations with a VOF method and carefully instrumented model-scale experiments are carried out to investigate extreme runup events around a ship-shaped body in transient wave groups. The numerical results are verified by comparing with experimental measurements at key locations, and are extended to explore global hydrodynamic behaviour across the full spatial domain. Thus, CFD-based numerical modelling is demonstrated to be an effective tool for exploring the physics of wave–structure interaction leading to the onset of greenwater on deck for a simplified FPSO geometry. Accurate predictions of the local wave field and the vessel motion are obtained.
Using the phase separation method, all the harmonic interactions up to (at least) the fourth are resolved. The extracted linear and high-frequency second harmonic components agree well with diffraction solutions up to second order. This indicates that the Stokes-like approximation is appropriate, and apart from the second-order difference terms, the subharmonic frequency components at each order in magnitude are relatively small even for the complex and violent interactions of a steep wave group with a fixed or a freely floating structure. The low-frequency second harmonic component is less well modelled, but this is attributed to ‘error waves’ in the CFD simulations.
This work is a generalization of the Stokes-wave perturbation expansion to wave runups around the vertical bow of a floating FPSO. When compared to Chen et al. (Reference Chen, Zang, Taylor, Sun, Morgan, Grice, Orszaghova and M Tello2018) and Chen et al. (Reference Chen, Stagonas, Santo, Buldakov, Simons, Taylor and Zang2019a), for which the phase separation method is applied to explore Stokes-type perturbation expansion to wave forces on a vertical cylinder, the scattering body in this work has a very different shape, and it is freely floating so is allowed to move under the action of waves. The methodology is not directly applicable to the analysis of wave runup at a flared bow, which may show a jet like feature.
The verified models/analysis are utilized to demonstrate that both low- and high-frequency second harmonic components would lead to wave runup at significantly higher levels than expected. Additionally, there are no large higher harmonics beyond second order. The nonlinearity in locally diffracted waves rather than vessel motion is the key mechanism for the excitation of nonlinear extreme runup. Ship motion must be included but a linear representation of this appears to be adequate. This is consistent with the observations in Stansberg & Berget (Reference Stansberg and Berget2009) and Schiller et al. (Reference Schiller, Pâkozdi, Stansberg, Yuba and e Silva2014).
The findings show that second-order diffraction theory may be accurate enough to justify development of an advanced engineering screening tool that would reduce statistical uncertainty in random sea states. At present, widely used screening tools for extreme runup, hence greenwater events, are generally based on linear wave theory (with/without an empirical correction). The existence of the bow flare or other changes in vessel geometry at and above the mean sea water level may play important roles and we leave investigation of these additional complications for future work. The possible application of the methodology to a wider range of wave conditions and vessel headings (of practical interest) will also be carried out and will be reported separately.
Supplementary movies
Supplementary movies are available at https://doi.org/10.1017/jfm.2020.1072.
Acknowledgements
The authors are grateful to Professor B. Teng at Dalian University of Technology for general discussions on the diffraction analysis and the access to the diffraction code. We thank Mr Y. Pei, H. Li, Q. Kuang, D. Mu and Y. Zhou and all the staff at Dalian University of Technology, China for their help and hospitality during the experiments.
Funding
This work was supported by the ARC Industrial Transformation Research Hub for Offshore Floating Facilities which is funded by the Australian Research Council, Woodside Energy, Shell, Bureau Veritas and Lloyds Register (grant no. IH140100012) and was supported by the National Natural Science Foundation of China (grant no. 52001053). D.Z.N. is grateful for support through Liaoning BaiQianWan Talents Program. The support in resources provided by the Pawsey Supercomputing Centre with funding from the Australian Government and the Government of Western Australia is acknowledged. Part of this work was supported by computational resources provided by the Australian Government through Gadi under the National Computational Merit Allocation Scheme.
Declaration of interests
The authors report no conflict of interest.
