## Abstract

A general phase-based harmonic separation method for the hydrodynamic loading on a fixed structure in water waves of moderate steepness is proposed. An existing method demonstrated in the experimental study described by Zang *et al.* (Zang *et al.* 2010 In *Proc. Third Int. Conf. on Appl. of Phys. Modelling to Port and Coastal Protection.* pp. 1–7.) achieves the separation of a total diffraction force into odd and even harmonics by controlling the phase of incident focused waves. Underlying this method is the assumption that the hydrodynamic force in focused waves possesses a Stokes-like structure. Under the same assumption, it is shown here how the harmonic separation method can be generalized, so that the first four sum harmonics can be separated by phase control and linear combinations of the resultant time-histories. The effectiveness of the method is demonstrated by comparisons of the Fourier transforms of the combined time-histories containing the harmonics of interest. The local wave elevations around the focus time are also visualized for the first three harmonics in order to reveal the local dynamics driving components within the wave force time-history.

## 1. Introduction

Offshore wind turbines are typically mounted on bottom-founded piles and are situated in areas of high winds where large waves can be expected to occur. The wave-induced loads on the foundations owing to extreme incident waves will comprise a linear harmonic component at around the incident spectral peak frequency and higher-order harmonics which are at multiples of the linear peak frequencies. These arise from the nonlinearity inherent in the waves and in the interaction. The higher-order components of the total wave force are of interest in order to predict the structural response. This is because the higher-order harmonics arising in the wave structure interaction can cause a large dynamic response of the piles owing to the excitation of a structural resonance in the lightly damped system. This resonant excitation can occur even when the high-frequency incident wave and local diffracted wave components are small relative to the total wave field. The excitation of sharp transient structural response by the higher-order harmonics (typically the component at the third-order sum frequency) is referred to as ‘ringing’ and has been demonstrated experimentally by among others [1,2]. Both studies involved vertical circular cylinders. A prior review by [3] discusses how ringing was observed to arise in steep wave interactions with fixed gravity-based structures and floating tension leg platforms at model and full scale. Thus, ringing has been a topic of interest within the field of offshore engineering for a significant period of time. Given the high-frequency nature of the ringing excitation, it is important to isolate and accurately compute the higher-order wave loading harmonics in order to fully capture the total forces and response of vertical cylinders in extreme waves. Rather than the phenomenon of ringing itself, it is the separation of higher-order force harmonics that is considered in detail here.

Extreme waves which arise during prolonged storms can be effectively modelled using focused wave groups and in particular the NewWave concept described by [4]. NewWave-type focused wave groups correspond to the average shape of the free surface in the vicinity of large crests occurring in a random sea state and are used throughout the analysis to follow. As shown by [5,6], the NewWave model captures quite accurately the linear contribution to large ocean waves. In this paper, we seek to understand the effects of nonlinearity in the interaction between focused waves and a surface-piercing cylinder on the wave-induced loads. The particular focused waves considered comprise a linear NewWave component with additional higher-order terms owing to wave–wave interactions. Experimental studies of the harmonic structure of the loads on a single surface-piercing column for incident focused wave groups of various degrees of nonlinearity have been conducted by [7,8]. In both cases, the challenge was to extract the higher-order harmonics without recourse to a Fourier transform of the force time-histories—an approach used in many of the preceding studies involving regular waves—which is not suitable for transient interactions.

By assuming the existence of a generalized Stokes-type harmonic series in both frequency and wave steepness for the free-surface elevation and horizontal wave loads for focused waves groups, Zang *et al.* [8] extracted the higher harmonics in the interaction. An alternative approach taken by [7] involved the use of a scaled-averaged wavelet energy analysis to isolate the harmonic force amplitudes. The dominant loads identified by [8] were found to be potential in form and, although the separation of drag and potential flow contributions to the total force was not possible, the effects of drag should be small based on Keulegan–Carpenter number estimates. The higher-order nonlinear contributions to the free-surface elevation and wave loads (also referred to as higher-order harmonics) were separated using the ‘phase-inversion’ method described by [6,9–11] for the wave kinematics at the focus point of the wave group. To extract the odd and even harmonic contributions to the time histories of kinematic or dynamic quantities (e.g. the free-surface elevation or total wave force) in the focused wave group interaction, it is necessary to conduct two tests involving incident wave groups with identical component amplitudes and frequencies but inverted phases. The individual harmonics can then be isolated by digital frequency filtering. The harmonic structure of the total horizontal force was also demonstrated in [8], even for violent breaking waves, and the amplitude and timescale of the higher-order harmonics were shown to correspond to the linear envelope/component raised to the appropriate power. Thus, the generalization of the Stokes-type expansion to wave groups was shown to be appropriate in describing the dynamics of the interaction.

Harmonic separation techniques can be applied to fully nonlinear time-histories arising from experimental tests or numerical simulations. The advantages of applying harmonic separation methods to numerical simulation data are that not only can the harmonics of the force time-histories be separated, but also the harmonics of the local flow field. Here, time-histories from fully nonlinear potential flow simulations of the experiments described by [8] involving focused wave diffraction by a vertical circular cylinder are analysed. The numerical scheme adopted for these fully nonlinear simulations was the mixed Eulerian–Lagrangian higher-order boundary element method (BEM) described by [12,13]. Other numerical solution approaches for the potential-flow diffraction problem in the nonlinear wave regime include a finite-element computation based on Hamilton's principle [14], and a nonlinear decomposition method to solve for the diffracted wave field assuming the incident wave is explicitly known as described by [15] and recently investigated in more detail in [16]. Viscous flow solvers have also been used in order to describe the physics of the interaction more comprehensively, including the hybrid finite volume–volume of fluid method adopted by [17] and a spectral wave explicit Navier–Stokes equations solver used by [18].

In this paper, we extend the existing phase-based separation method (i.e. the ‘phase-inversion’ method), so that the first four harmonics of the force or elevation can be separated by the application of simple linear combinations of the relevant time-histories with minimal post-processing. This requires the data from four diffraction interactions involving incident-focused waves that are each exactly *π*/2 out of phase. Therefore, fully nonlinear potential flow simulations of two focused-wave interactions not considered experimentally are also required. The advantages of the generalized four-phase method compared with the existing two-phase method are demonstrated by the comparison of a higher harmonic signal obtained from each method. A comparison of the four-phase harmonic separation results for numerical time-histories to the results of the phase-inversion method applied to experimental force time-histories is made in order to validate the numerical approach and to confirm convergence. The four-phase harmonic separation method is also applied to the elevation of the free-surface surrounding the body and to the time-history of the wave run-up on the cylinder. The evolution of the local flow field around the body during the interaction at first, second and third harmonic is presented and discussed.

The structure of the paper is as follows. In §2, the extension of the Stokes wave expansion from regular waves to wave groups is discussed. Thereafter, in §3, the extension of the phase-inversion method, which separates odd and even harmonics, to a more general phase-based separation method that isolates the first four harmonics is presented. The application of this generalized phase-based separation method to fully nonlinear simulations of an experiment is described in §4 with comparisons with these experimental results also provided. The behaviour of the first three harmonics of the local flow field during the diffraction interaction is illustrated and discussed in §5. Finally, a conclusion regarding the advantages of the generalized phase-based separation is presented in §6.

## 2. Stokes wave expansion generalized to wave groups

The interaction under consideration involves the diffraction of a unidirectional incident focused wave by a surface-piercing bottom-founded cylinder. It is assumed that the viscous effects are negligible relative to the potential flow effects, and so the equations for an irrotational flow in an inviscid, incompressible fluid are used to describe the diffraction. The implications of this assumption are considered in §4*a*,*c*. According to the classic Stokes perturbation expansion in the wave slope parameter *kA*, a regular incident wave of linear wave amplitude *A*, frequency *ω* and wavenumber *k* will induce a horizontal hydrodynamic load of the form
*A*, where the coefficients *f*_{mn} represent wave-to-force transfer functions, and *φ*=*ωt*+*φ*_{0} is the phase of the linear component of the incident wave with a prescribed phase shift *φ*_{0}. In order to keep this force expression as simple and compact as possible, the term *Re*{*f*_{mn} *e*^{inφ+iψmn}}, where *ψ*_{mn} is the phase shift of the force relative to incident wave at the given order and harmonic. The structure of the terms multiplying the coefficients *f*_{mm}, where the power of the amplitude term is the same as the frequency harmonic, and the terms multiplying coefficients such as *f*_{m(m−2)}, where the power of the amplitude exceeds that of the frequency harmonic by 2, is noteworthy. In the notation of Stokes' expansions, these coefficients correspond to sum and difference frequency components, respectively. For ringing, we are mostly concerned with the sum harmonics.

While these individual harmonic components are easy to extract for regular wave forcing, this is much more difficult for broad-banded waves trains (either random or wave groups) where the simple Stokes' terms are replaced by summations of products of linear terms, and each net higher harmonic contributes across an increasingly broad range of frequencies. This broadening of the higher harmonic spectral peaks is referred to as spectral ‘smearing’. For example, the Stokes' perturbation expansion to second order in *kA* for the free-surface elevation in an irregular sea can be expressed as
*A*_{n} is the amplitude of the free-wave component of frequency *ω*_{n}, wavenumber *k*_{n} and phase *θ*_{n}=*k*_{n}*x*−*ω*_{n}*t*. The coefficients *B*^{+}_{mn} and *B*^{−}_{mn} account for the super- and subharmonic contributions at second-order arising from wave–wave interactions and are given by [19] for constant finite depth and are referred to as bound waves. Analytical third-order solutions for irregular waves in deep water have also been derived by [20]. The wave loading on a cylinder generalized to an irregular sea can be written in a similar manner to equation (2.2) by replacing the coefficients *B*^{+}_{mn} and *B*^{−}_{mn} with wave-to-force quadratic transfer functions (with extensions to higher order).

For a focused wave group, the phases of each free-wave component are chosen, so that the crests of all linear components coincide at the focus event (*x*_{0},*t*_{0}), and significant nonlinearity will be confined to the neighbourhood of the focus location around the focus time. If a cylinder is placed at the focus location, then the wave loading will contain only significant nonlinearity in a short time-interval around the focus time. We propose that the regular wave Stokes equation can also be used to approximate the (focused) wave group case where the simple Stokes terms now represent summations of products (of linear terms). However, the convergence of the Stokes-type expansion for wave groups is not guaranteed. Nevertheless, we seek to extract the higher-order components as far as fourth order using phase manipulations to separate each component.

To demonstrate how the regular-wave Stokes' expansion might be generalized to wave groups, it is useful to assume the energy spectrum of the NewWave is such that the narrow-banded approximation can be adopted. Therefore, around the focus point, the signal envelope can be regarded as a slowly varying function of time *A*(*t*) (relative to the ‘carrier wave’ oscillations at the peak spectral frequency *ω*_{0}). In this case, the free-surface elevation to second order at the focus location around the focus time can be simplified to
*φ*_{0} is the phase of the peak frequency component and, *B*^{+} and *B*^{−} are the second-order sum and difference coefficients at the peak frequency *ω*_{0}. The narrow-banded approximation reduces the complexity of the equations significantly by approximating the summations and double summations of equation (2.2) with the product of a wave group envelope *A*(*t*) and the oscillations within the envelope *A* with time-dependent amplitude *A*(*t*) describing the envelope of the oscillations of peak frequency *ω*_{0}. For each higher harmonic, the contribution to the total force is expected to decrease in both amplitude and temporal extent owing to the *A*^{n}(*t*) term. In mathematical terms, the Stokes' expansion for a broad-banded wave group ultimately will not converge, because the contributions at higher harmonics eventually cease to decrease with increasing order. However, typical energy spectra used in experimental analyses (JONSWAP, Bretschneider) are sufficiently narrow-banded to assume the Stokes structure holds up to fifth order at least as shown by [8].

## 3. Generalization of the harmonic decomposition method

### (a) Harmonic decomposition via the phase-inversion method

The phase-inversion method applied to wave kinematics described by [10] can be used to separate the odd and even harmonics in the diffraction force time-history. In brief, it requires two time-histories for a given hydrodynamic quantity from two interactions involving incident waves generated by linear paddle signals that are 180° out of phase. By consideration of equation (2.1), it is clear that the phase shift *φ*_{0} is chosen, so that it generates a wave with a crest at the focus, and thus the inverted signal results in a trough at the focus. The even and odd harmonics of the fully nonlinear force time-history are then separated by frequency filtering. However, there are strict limits as to what can be achieved by filtering in frequency for higher-order harmonics owing to the ‘smearing’ of the spectrum as a result of wave–wave interactions. In particular, at higher orders, overlap between the *n*th and (*n*+2)th harmonics can occur over a range of frequencies. It should be noted that the Stokes expansion for the force is expressed in regular wave form but is representative of the generalized expansion for a wave group.

### (b) Generalization of the phase-inversion method

In order to obtain more effective separation of the harmonics, it is necessary to generate more force signals corresponding to further phase shifts of the linear paddle signal. Consider four diffraction interactions involving wave groups generated by the same paddle signal, except the phase of each Fourier component is shifted 0°, 90°, 180° and 270° leading to four force time-histories *F*_{0},*F*_{90},*F*_{180},*F*_{270}. We then seek to extract the individual Stokes-type components through linear combinations of these time-histories and the Hilbert transforms (i.e. the harmonic conjugate of the signal denoted by the superscript *H*) of the 90° and 270° force time-histories, rather than depending on frequency separation of the signals. The linear harmonic and the first three superharmonics (to fourth order) are of most interest and are sought in isolation from each other—the second-order subharmonic will occur in the same expression as one of the superharmonics. The required simple linear combinations of the four phase-shifted runs are as follows
*n*th and *n*+2th harmonics in the phase inversion method so that separating these contributions by digital filtering is trivial.) Thus, running four phase-shifted cases yields the dominant harmonic components up to the fourth sum harmonic. We note in passing the relative importance of the difference terms, such as *f*_{31} compared with the *f*_{11} term. These terms have the same frequency content, but a different (higher-order) dependency on the wave amplitude. In general, all such difference terms are likely to be negligible for weakly nonlinear waves, except for the second-order difference long-wave *f*_{20} term. Therefore, this harmonic separation method will only yield the linear harmonic and the first three superharmonics (with trivial filtering necessary only for the fourth-order superharmonic), if the waves are at most of moderate steepness. We stress that other phase-shifted combinations of wave groups would allow separation of the harmonics. The experimental work of Hann *et al.* [21] makes use of 12 versions of the same wave group, each shifted by 30°. It is clear that the four linear combinations that we use are the minimum number required to give clean separation of the first four harmonics.

## 4. Experimental and simulation results

The harmonic decomposition method generalized to four phase-shifted incident wave groups is demonstrated using time-histories obtained from potential flow simulations of a set of focused wave experiments. In particular, a fully nonlinear higher-order BEM potential flow model (described by [13]) is used to simulate the experimental tests involving isolated compact wave groups (with phase shifts of 0° and 180° from the crest focused case) incident on a single surface-piercing column reported in [8]. In addition to the crest and trough focused wave group experimental tests, we also simulated interactions involving wave groups with Fourier components shifted by 90° and 270° relative to the crest focused wave groups. A comparison with the harmonic separation obtained from the analysis of the experimental results (presented by [8,22]) is also provided.

### (a) Experimental tests: simulation case study

The numerical simulations use the experimental specifications for the smallest focused wave group considered, which are as follows: the cylinder has a diameter of 0.25 m, the distance from the wavemaker to the focus point at the front stagnation point of the cylinder is 7.8 m, and the water has a depth of 0.505 m. The particular test wave that is considered here is a unidirectional NewWave-type wave group with an underlying JONSWAP (*γ*=3.3) spectrum of peak frequency *f*_{p}=0.61 Hz. The linear component of the focused wave was prescribed to have an amplitude *A*=0.06 m corresponding to a steepness of *k*_{p}*A*≃0.1 (the smallest steepness in the programme of experiments), where *k*_{p} is the wavenumber corresponding to *f*_{p}. In the experiments, the wavemaker comprises a segmented piston paddle array, and a linear signal is used to generate the waves. In the simulations, a piston motion of a vertical wall is implemented to model the experimental wave generation. Good agreement between experimentally and numerically generated focused waves was observed. Implementing phase shifts in the numerical wavemaker signal and hence obtaining the phase-shifted force time-histories was straightforward.

### (b) Phase-based separation results

Four distinct simulations of the diffraction interaction were conducted each involving an incident wave phase-shifted by 90° relative to the preceding simulation. A Hilbert transform was applied to two of the time-histories, so that all the necessary data for obtaining the individual superharmonics from linear combinations (3.3)–(3.6) were available. A simple high-pass filter was applied to isolate the fourth-order superharmonic from the second-order subharmonic in equation (3.6).

To illustrate the success of the harmonic separation, it is useful to apply a Fourier transform to the four composite signals (3.3)–(3.6) obtained from the four linear combinations of the time histories and Hilbert transforms. Given the large difference in the amplitudes of the harmonics (from approx. 50 Newtons at first order to 1 Newton at fourth order), it is necessary for illustrative purposes to use a log-scaling of the Fourier transform amplitudes. The semi-log plot of the Fourier transform of the four composite force time-histories is shown in figure 1*a* (the FTs are relatively noisy, because no smoothing is applied to the data). The degree of separation between the harmonics is excellent. The semi-log plot of the Fourier transform of the odd and even harmonics separated using the two-phase approach (see equations (3.1) and (3.2)) is shown for comparison in figure 1*b*.

Although the signal is relatively noisy below an amplitude of 10^{−2}, the harmonic structure is visible up to fifth-order sum frequency as is evident from the secondary peak in the signal dominated by the linear harmonic at *f*=5*f*_{p}. This corresponds to the *O*(*a*^{5}) term in the first expression in equation (3.3) which can be shown to be the fifth harmonic *b* is relatively straightforward, because the overlap between harmonics is minimal—this is evident from a comparison of the amplitude spectrum for the even harmonic to those of the second-order sum harmonic and second-order difference and fourth-order sum harmonics in figure 1*a*. However, isolating the third-order harmonic is not as straightforward, because the amplitudes of the high-frequency ‘tail’ contributions of the first harmonic spectrum are of the same order of magnitude as the third harmonic peak contribution in figure 1*a* and occur in the same frequency range.

To isolate the third harmonic from the Fourier transform of the odd harmonic signal, it is necessary to apply a suitable filter. In §7.3 of [22], two possible filters for the *n*th sum harmonic are proposed: one of width *f*_{p} centred on *f*=*nf*_{p} and the other of width 2*f*_{p} centred on *f*=*nf*_{p}. Here, both narrow and wide filters are applied to the odd harmonic amplitude spectrum, and the corresponding time-signals are cross-compared with the third-order signal from the four-phase harmonic separation method. A representation of the digital filter applied to the amplitude spectrum of the odd harmonic signal (3.1) extracted using the two-phase approach is shown in figure 2*a*—the amplitude spectrum of the third harmonic linear combination (3.5) from the four-phase method is also plotted in figure 2. The associated time-histories are compared in figure 2*b*. From this comparison, it is evident that the narrow filter does not ‘pass’ enough of the third-order contribution, in particular the tail of the third-order peak, so that the amplitude of the force oscillation is smaller than for the four-phase case. By contrast, the wide filter includes contributions from the odd harmonic spectrum over an excessive frequency range including some of the tail of the first-order peak. The resultant signal is too large relative to the four-phase third harmonic signal. The spectral overlap of the first and third harmonic contributions is unavoidable, and so the third harmonic cannot be isolated from the odd harmonic spectrum obtained using the two-phase separation method. Of the two filters, the narrow version yields a more accurate representation of the third-order contribution and the same approach is used for all the higher-order harmonics in the phase inversion method.

A comparison of the first four sum harmonics and the second-order difference harmonic as obtained from the numerical simulations using the four-phase combination and the two-phase combination (narrow filter) methods is shown in figure 3. There is a significant difference in the shape of the wave packets obtained using the two- and four-phase approaches for the second-order difference, third and fourth harmonics. In the case of the third and fourth harmonics, the amplitude of the oscillations at the focus time is smaller for the two-phase approach, and the oscillations prior to the focus time have a slightly larger amplitude compared with the four-phase approach. For the first two harmonics, the differences between the time-histories extracted using the different phase separation approaches are minimal.

It is useful to note that if the wide filter had been used, then the third- and fourth-order results would be larger for the two-phase approach than for the four-phase approach. The unambiguous time-histories obtained from the four-phase harmonic separation method contrast to the filter-dependent time-histories obtained from the two-phase method. Furthermore, assuming (as often happens when ringing arises in diffraction problems) the third-order harmonic force excites the structural resonance, then the accurate isolation of this harmonic contribution will be of critical importance. The four-phase harmonic separation method returns an unambiguous third-order sum harmonic time-history, whereas the third-order signal depends on the range of the spectral filter in the two-phase method. Therefore, the four-phase harmonic separation method can be considered a significant improvement on the two-phase method when investigating the excitation of high-frequency ringing loads on structures.

### (c) Computational mesh specifications and preliminary validation

The amplitude spectra of the force harmonics shown in figure 1 indicate that the numerical simulation also resolves the fully nonlinear interaction up to fourth order with possible fifth-order contributions. The properties of the mesh were specified in order to resolve the fourth-order waves while avoiding prohibitive simulation computational times. A proper convergence investigation would require a simulation involving a finer mesh. However, the results presented correspond to approximately 6300 quadratic boundary elements and 15 800 nodes in the numerical model. The computational time for the simulations was approximately 800 CPU hours, and a significant increase in the number of elements/nodes would have resulted in excessively long computational times.

Nevertheless, it is useful to describe the main properties of the computational mesh implemented in the fully nonlinear simulations. The unstructured free-surface mesh of triangular elements implemented in the BEM as part of the fully nonlinear simulation (see [13] for further details) are of side-length 0.089 (where one unit of length equals a metre in the numerical simulation). The wavelength of the peak frequency component of the incident focused wave is 3.19 m, and so there are approximately 36 free-surface elements per first-order wavelength. For the fourth-order bound waves of wavelength 0.80 m, there are approximately eight elements per wavelength and for the fourth-order free-waves (0.26 m) only three elements per wavelength. The boundary data (position coordinate, velocity potential) vary quadratically on the elements and on each triangular element there are six nodes. Such a higher-order BEM does not require a large density of elements to resolve a given wave component. Furthermore, immediately around the body where free-waves will be generated the free-surface is finer with mesh elements at the cylinder of length 0.033 corresponding to eight elements per fourth-order free-wavelength. This mesh was considered sufficiently fine to model the nonlinear interaction up to fourth order. The time-step was chosen to be 0.01 s corresponding to 1/160th of the peak period.

A preliminary validation of the simulation results was possible by comparing the results of harmonic decomposition of the numerical time-histories to the diffraction load harmonics measured in the experimental tests (decomposed as described by [8]). This comparison of the numerical results to the experimental results is also shown in figure 3. Apart from the third-order sum harmonic, the agreement between the simulations and experimental tests is very good. There are a number of possible reasons for the large discrepancy between the numerical and experimental results at third order. First, viscous effects such as drag- and wave-induced vortex shedding may play a role in the higher-order harmonics. If the representative time scale of the interaction is considered to be the peak period, then the Keulegan–Carpenter parameter for the interaction is *KC*≃1.5. Therefore, the drag contribution will be significantly smaller than the potential flow contribution to the total load, but may be large enough to affect the higher harmonic forces such as the third harmonic (large-scale wave-induced vortex shedding is unlikely to occur for this *KC* value). Another possible reason for the discrepancy is the use of the two-phase harmonic separation method for processing the experimental data—the filtering of the amplitude spectrum may not have completely isolated the third-order contribution from the first-order contribution. Although the amplitude of the fourth-order time-histories do not agree exactly just after the focus time *t*=0, the phase of the oscillations aligns almost exactly. Therefore, it is concluded that the interaction is simulated accurately to fourth order.

It is beneficial to explore the drag contribution to the experimental force in more detail in order to understand the applicability of the model to experimental diffraction forces. The Morison drag term in a sinusoidal flow field can be expanded into harmonic components as *O*(*a*^{2}) is likely to be much smaller than the corresponding potential flow force harmonic which is of order *O*(*a*); however, the third harmonic potential flow force is *O*(*a*^{3}) and so smaller than

## 5. Numerical solution: analysis and discussion

The four-phase harmonic separation method summarized by equations (3.3)–(3.6) allows the properties of the solution at the individual-separated harmonics—such as the diffraction force time-history and the wave field surrounding the body—to be investigated in a straightforward manner. To extract the first four harmonics from the fully nonlinear diffraction force time-history (given four fully nonlinear time-histories for focused waves 90° out of phase), it is necessary to apply the Hilbert transform to two of the nonlinear time-histories and to filter equation (3.3) as previously described. This procedure requires minimal computational effort. However, the separation of the harmonics of the nonlinear free-surface throughout the computational domain is much more computationally intensive as it necessitates taking the Hilbert transform of the elevation of each node on the free-surface over a specified duration of time. Given that the free-surface mesh changes in time, this requires that the elevation of the free-surface nodes at each output time must be interpolated onto a fixed mesh (e.g. the mesh at the first time-step). Thereafter, the elevation time-history of each of the nodes on this fixed mesh must be processed using the Hilbert transform. Once these steps have been taken, the harmonic separation equations may be applied to obtain the separated harmonics at each instant at which the free-surface elevation is outputted.

### (a) Harmonic force time-histories

The harmonic separation of the diffraction force time-history has already been described and discussed (albeit briefly) in §4*b*. However, to put the higher-order harmonic force time-histories into context, it is useful to compare the harmonics with the powers of the linear harmonic envelope. In [8], it is noted that the envelope of the first five higher-order harmonic contributions in time match the envelope of the linear force component raised to the appropriate power and scaled according to the magnitude of the appropriate harmonic. The envelope of the linear force component *F*_{1} is obtained from *F*^{H}_{1} is the Hilbert transform of the first-order force time-history. Therefore, in obtaining the estimated wave envelope for the third-order harmonic, the quantity

By contrast, the third harmonic force time-history displays a distinct behaviour of its own. The oscillations in the force time history do not follow the linear envelope raised to the power of three, rather there is a delay relative to the envelope prediction before significant oscillations occur. Furthermore, the third harmonic oscillations persist with an amplitude close to the peak amplitude for a slightly longer time compared with the envelope prediction. This indicates that some particular interaction mechanisms arise at third order that do not arise at second or fourth order.

To further examine the properties of the third harmonic force time-history, the contributions from the third harmonic amplitude spectrum immediately around that spectrum's peak (*f*≃3*f*_{p}) are separated from the ‘high-frequency’ tail (*f*>3.5*f*_{p}) contributions. These contributions are then compared with the high-frequency tail component of the second harmonic and the peak fourth harmonic component, respectively. The second, third and fourth harmonic amplitude spectra are all cleanly separated; we are interested in comparing the behaviour of the components of the amplitude spectra on the same frequency ranges. These contributions to the isolated harmonics are extracted using a simple frequency filter. It must be emphasized that this filter is required only to examine the properties of the third harmonic—it has not been used in producing the spectra shown in figure 1*a* or for the numerical results in figure 3. The filter is designed to fully pass the spectral contributions in the frequency interval (2.75*f*_{p},3.25*f*_{p}) and to linearly ramp the contributions to zero from the 2.75*f*_{p} to 2.5*f*_{p} and from 3.25*f*_{p} to 3.5*f*_{p}. It is clear from figure 5*a* that the high-frequency ‘tail’ of the second harmonic and the peak third harmonic contributions correspond to forces of a quite different time-dependence. The maximum force of this component of the third harmonic occurs approximately 1 s later than the high-frequency component of the second harmonic force—this accounts for the delayed maximum observed in figure 4*b*. It is not clear why the occurrence of the maximum amplitude is delayed. However, the harmonic separation method coupled with some frequency filtering has allowed a more detailed examination of the third harmonic. A similar comparison of the forces corresponding with the high-frequency tail component of the third harmonic and the peak frequency component of the fourth harmonic for the frequency range (3.75*f*_{p},4.25*f*_{p}) (with linearly ramped spectral contributions to either side of this interval included as before) is shown in figure 5*b*. In this case, the force time-histories are very similar—the maxima occur at approximately the same time—but the oscillations of the high-frequency spectral contribution to the third harmonic force begin slightly later than the oscillations corresponding to the overlapping fourth harmonic spectral contribution.

Therefore, the force time-history at the third harmonic displays different properties to those at the second and fourth harmonic. To identify possible differences between the interactions at each harmonic, it is useful consider the evolution of the wave field around the cylinder at each harmonic. To complement this decomposition of the free-surface elevation around the body into individual harmonics, the first three harmonics of the nonlinear run-up on the cylinder are also extracted and presented.

### (b) Scattered wave field

To investigate the diffraction processes that yield the particular force time-histories described above, it is useful to plot the total and diffracted wave field during the interaction. To obtain the diffracted wave-field pattern from the numerical simulations it is necessary to simulate both the diffraction problem with the cylinder present and the incident wave-propagation problem in the absence of the cylinder. The incident wave field is obtained using simulations involving a narrow numerical wave tank (comprising a free-surface two elements wide) thus avoiding long computational times. The free-surface elevation on the centre line is used to determine the two-dimensional wave field. Given that the free-surface meshes for the propagation and diffraction are different, the incident wave field must first be interpolated onto the (initial) free-surface mesh from the diffraction problem before subtracting the incident field from the total field to yield the diffracted wave field. This process is repeated for the first three harmonics only. To obtain the fourth-order-diffracted wave field, it is necessary to filter the free-surface elevation time-history obtained from equation (3.6) to extract the sum frequency term for every node on the free surface. Such a process would be quite time-consuming and computationally demanding and, so only the second and third harmonic wave fields were considered.

The features of the linear-, second- and third-order harmonic wave fields are highlighted by plotting sequences of contour plots for the free-surface elevations around the cylinder as shown in figures 6, 7 and 8, respectively. In figures 6–8, the upper half of the plots illustrates the contours of the total wave field and the lower half illustrates the contours of the diffracted wave field. Although the computational domain has a half-width of 1.5 m and length approximately 10 times the half-width, only a reduced square subdomain in the immediate vicinity of the cylinder is plotted for each harmonic (with the subdomain extent smaller for the higher harmonics) in order to emphasize the local wave field structure. The free-surface is plotted at different instants around the focus time *t*=0 s with time-increment being larger for the linear harmonic free-surface plots (Δ*t*=0.20 s) than for the corresponding second and third harmonic plots. To put these time-steps in context, it is reiterated that the peak period for the incident focused waves is 1.64 s.

A cursory inspection of the contour plots reveals that the diffracted components of the free-surface elevation at second and third order are significantly larger (relative to the incident component of the free-surface) and more complex in form than for the linear harmonic. For the linear component of the interaction, the peak incident wavelength λ_{1}≃3.2 m is much larger than the cylinder dimensions (diameter *D*=0.25 m), so that the diffracted wave field is relatively small and exhibits no remarkable properties beyond alternating positive and negative diffracted wave run-up at the front and rear face of the cylinder. The dominance of the incident wave implies the first harmonic can be well approximated by the Morison inertia force. However, the peak second-order-bound incident wavelength is half that of the first-order wave λ_{2}≃1.6 m, so the diffracted wave field and incident wave field are of similar order, as is clear from figure 7. The same argument applies for the third harmonic component of the free-surface. In fact, for the third harmonic, the local diffracted wave field is significantly larger than the incident wave component as is clear from figure 8 (it may be recalled that the upper figures show the sum of incident and diffracted waves). The local structure of the free-surface immediately surrounding the cylinder for the second and third harmonics features distinct diffraction patterns which show the generation of free-waves which propagate laterally outwards from the cylinder. The structure of the diffracted wave is quite different in each case.

In the contour plots for the second (sum) harmonic component of the free-surface elevation shown in figure 7, there is a clear diffracted wave component detaching from the cylinder and radiating outwards in the approximate direction *θ*=±*π*/2 radians (*θ*=0 corresponds to the direction of incident wave propagation). A large diffracted wave is also reflected back towards the wavemaker and spreads into a large arc some of which is evident in the contour plot at *t*=0.2 s and *t*=0.3 s. This wave arises owing to the large run-up and local free-surface distortion that occurs at the front of the cylinder. A small wave is also radiated from the downstream stagnation point but is not as significant as either the radiated wave travelling upstream or the radiated wave detaching from the side of the cylinder.

A comparison of figure 8 with figure 6 shows how the local third harmonic diffracted wave field is much more significant than the third harmonic incident wave field in contrast to the dominance of the incident wave field at first order. (Nonetheless, in absolute terms, the third-order-diffracted wave amplitude is smaller than that of the first-order-diffracted wave amplitude.) Furthermore, the free-surface variations at third order are of much more spatially compact than at first order which is why the contour plots are shown over a smaller region than the contour plots of the former. The most striking aspect of the third-order-diffracted wave field is the lateral radiation of two different wave disturbances from each side of the cylinder—as opposed to the single wave disturbance radiated at second order. On the downstream side of the cylinder, a wave disturbance detaches from the cylinder and radiates outwards at an angle approximately equal to *θ*=±*π*/3. On the upstream side of the cylinder wave diffraction occurs through the radiation of a free wave in the direction *θ*≃2*π*/3. Wave scattering from the upstream stagnation point is also visible although the wave amplitude is smaller than in the second harmonic case. The evolution of the upstream diffracted wave field and, in particular, the final three snapshots in figure 8 bear some resemblance to the scattered waves reported by [23]. The high-frequency-scattered waves documented in that paper were observed for cylinder diameter-to-wavelength ratios similar to the case considered here. However, rather than focused waves, steep regular incident waves were examined in more detail in that paper. The generation of the free waves at the two lateral points (as opposed to the front and rear of the cylinder) on the cylinder circumference may have an influence on the third harmonic force time-history and its divergence from the first-order envelope prediction. The distance between a successive positive and negative radiated disturbance is a little greater than 0.2 m and this half-wavelength is approximately equal to the peak wavelength of the third harmonic free-waves computed to be 0.46 m. It should be noted that this wave diffraction pattern persists well beyond the time-interval illustrated in figure 8 and is visible past *t*=1.0 s.

### (c) Wave run-up on cylinder

The fully nonlinear wave run-up on the cylinder is also decomposed into its constituent harmonics. Rather than considering the total run-up, however, the run-up owing to diffraction only is considered. The time-history of the first three harmonics of the wave run-up on the cylinder is represented through the (*t*,*θ*) surfaces illustrated in figures 9*a*,*b* and 10. The first harmonic of the diffracted wave run-up on the cylinder is shown in figure 9*a,* and the presence of the cylinder induces an alternating magnification and reduction of the free-surface elevation at the upstream and downstream face of the cylinder, respectively. Thus, for *t*≃−0.4 s, the diffracted wave yields a significant increase in the run-up on the cylinder around the front stagnation point (*θ*=180°) while simultaneously yielding significant decrease in the run-up around the rear stagnation point (*θ*=0°). (The increase and decrease is relative to the incident wave free-surface elevation at the cylinder.) No significant diffracted wave run-up is observed at the cylinder ‘shoulder’, i.e. at *θ*=90°. This behaviour is consistent with the flow-field contour plots in figure 6 and, more importantly, with the nonlinear wave slender body theory described by [24] generalized for irregular waves by [25]. In this case, the first harmonic component of the diffracted wave elevation at the cylinder is given by
*A* is the slowly varying linear focused wave envelope *A*(*t*), and *ω*_{0} is the peak angular frequency with associated peak wavenumber *k*_{0}. The free-surface elevation at the cylinder as computed by this slender body approximation is shown in figure 9*c* and it compares very well with the first harmonic of the nonlinear potential flow results.

A similar comparison is made for the second harmonic component of the diffracted run-up. First, it should be noted that the presence of a large diffracted wave disturbance at the shoulder of the cylinder in the fully nonlinear model computations of wave run-up (figure 9*b*) is consistent with the second harmonic local flow field shown in figure 7. Furthermore, the perturbation analysis of [24] (‘FNV’ theory) can be used to obtain second harmonic free-surface elevation to *O*(*A*^{2}):
*ζ*_{2} is evaluated on the cylinder *r*=*a*, and the narrow-banded approximation has been applied in the same manner as for equation (5.1). The FNV theory prediction for the second harmonic component of the diffracted wave run-up shown captures, qualitatively at least, much of the run-up behaviour present in the fully nonlinear computational results as can be seen from a comparison of figure 9*d* to *b*. In particular, the large free-surface oscillations at *θ*=90° are present in the FNV theory in addition to some small oscillations at the downstream stagnation point. However, the significant run-up occurring at the upstream face of the cylinder is not present in the FNV run-up prediction. Furthermore, the FNV theory yields a significantly smaller amplitude for the second harmonic run-up as can be seen from the elevation scales in each case. These differences between the theory and the fully nonlinear results arise, because the cylinder does not satisfy the slender body assumption *ka*≪1, assumed to hold for all harmonics, at the second harmonic. For the first harmonic, the wavenumber of the peak incident wave component corresponds to a slenderness parameter *ka*≃0.25, so the cylinder is approximately slender. However, for the second harmonic, the bound (peak) incident wavelength is such that *ka*≃0.5 and thus the cylinder dimensions are effectively larger, and the slender body theory is less applicable. Therefore, it might be expected that a larger wave run-up would occur at the front of the cylinder in the fully nonlinear solution than that predicted by the FNV theory.

The behaviour of the third-order harmonic component of the diffracted wave run-up, shown in figure 10*b*, is not compared with the FNV theory for two reasons. First, an explicit third-order perturbation expression is not presented by [24] for the free-surface elevation. Second, although it is possible to obtain such an expression, the assumptions behind the theory will be violated at the third harmonic and even qualitative agreement with the full solution is unlikely. Nevertheless, the diffracted wave run-up at third harmonic exhibits a behaviour consistent with the corresponding local diffracted wave field plots in figure 8. That is, the run-up is minimal at the shoulder and maximal at approximately *θ*=*π*/3 and *θ*=2*π*/3. These diffracted wave disturbances are also clear in figure 8 where radiation of these disturbances occurs periodically. An important distinction between the second- and third-order-diffracted wave run-up is that at third order the run-up pattern is not symmetric around the focus time—it is significant for *t*≥−0.75 s until *t*≤1.5 s, whereas the second harmonic component is approximately symmetric in time. Such behaviour confirms the observations of the second and third harmonic force time-histories. The temporal structure of the third harmonic wave elevation and force appears quite different than that of the second harmonic (and fourth harmonic for the force as is evident in figure 4). This suggests that the local dynamics (which may give rise to the free wave structure around the cylinder circumference) are different at the third harmonic. The occurrence of local dynamics particular to the third harmonic may also explain the apparent dual nature of the third harmonic force, i.e. the different properties of the two third harmonic force components illustrated in figure 5. A different treatment of the third-order contributions to the force on a slender cylinder in nonlinear waves compared with the first and second order is required in the FNV analysis. Thus, the introduction of the third-order nonlinear potential *Ψ* in this analysis of the solution in the inner domain suggests that new sources of nonlinear dynamics arise beyond second order. Nevertheless, the fourth harmonic force time-history appears consistent with that of the second harmonic.

## 6. Conclusion

A generalization of the phase-based separation method for fully nonlinear free-surface elevation or hydrodynamic forces, used previously by [10] among others, has been presented here. The effectiveness of this method, which requires four time-histories obtained during interactions involving four focused waves 90° out of phase, has been demonstrated for the case of a JONSWAP NewWave incident on a single cylinder. The particular case considered used a fully nonlinear potential flow code to generate the time-histories; however, the method is applicable to both numerical and experimental investigations. In the latter case, the method requires the wave steepnesses to be at most moderate and viscous drag effects to be small. Separation of the first four superharmonics and the second-order subharmonic can be achieved with simple linear combinations of the time-histories. In experimental studies where focused wave group simulations are extremely short in duration, minimal extra effort is needed to generate four time-histories rather than one time history for a given wave. Therefore, this method is considered to be very useful for experimental investigations involving focused waves based on typical energy spectra (e.g. JONSWAP, Bretschneider) where higher-order nonlinear contributions are of interest. The method is not specific to the wave packet considered here and does not require the spectrum to be narrow-banded.

A particular advantage to using the generalized four-phase harmonic separation method is that only a trivial frequency filter is necessary to isolate the second-order difference harmonic and fourth-order sum harmonic. The first-, second- and third-order sum harmonics are obtained directly from the linear combinations of the four force time-histories (each *π*/2 out of phase) and the two Hilbert transform time histories. By contrast, careful filtering of the odd and even harmonic force time-histories is necessary in the two-phase harmonic separation method where spectral ‘leakage’ between harmonics can occur. This is particularly problematic when isolating the third harmonic in the two-phase approach owing to the disparate magnitudes of the first and third harmonics wherein the high-frequency tail of the first harmonic spectrum overlaps and has a similar order of magnitude to the third harmonic spectral peak. For investigations into the phenomenon of ringing the third-order sum frequency can often be the most important higher-order term and thus isolating the contribution effectively is important.

The structure of the wave fields and the wave run-up at each of the first three harmonics (linear, second-order sum, third-order sum) was also investigated through contour plots of the corresponding free-surface elevations obtained from the linear combinations in the four-phase harmonic separation method. Free wave radiation from the cylinder in one predominant direction was observed in the second harmonic of the diffracted wave field. Similar behaviour was noted by [26] in a second-order analysis of wave diffraction by a cylinder. The fully nonlinear simulation results allow us also to analyse the third harmonic of the diffracted wave field. For this third harmonic, two different wave disturbances were observed to radiate from points on the circumference of the cylinder located at *θ*=*π*/4 (leading to lateral wave radiation in the downstream direction) and *θ*=3*π*/4 (leading to lateral wave radiation in the upstream direction). Furthermore, qualitative agreement was observed between the fully nonlinear computational results for the diffracted wave run-up and run-up predictions from an analytical nonlinear wave slender body theory.

## Acknowledgements

This work was undertaken as part of the PerAWaT project, which was commissioned and supported by the Energy Technologies Institute. Numerical results were obtained with the use of the Oxford Supercomputing Centre facilities.

- Received December 20, 2013.
- Accepted May 20, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.