## Abstract

A new stochastic representation of a seastate is developed based on the Karhunen–Loeve spectral decomposition of stochastic signals and the use of Slepian prolate spheroidal wave functions with a tunable bandwidth parameter. The new representation allows the description of stochastic ocean waves in terms of a few independent sources of uncertainty when the traditional representation of a seastate in terms of Fourier series requires an order of magnitude more independent components. The new representation leads to parsimonious stochastic models of the ambient wave kinematics and of the nonlinear loads and responses of ships and offshore platforms. The use of the new representation is discussed for the derivation of critical wave episodes, the derivation of up-crossing rates of nonlinear loads and responses and the joint stochastic representation of correlated wave and wind profiles for use in the design of fixed or floating offshore wind turbines. The forecasting is also discussed of wave elevation records and vessel responses for use in energy yield enhancement of compliant floating wind turbines.

## 1. Introduction

The design of ships and offshore structures that are expected to withstand severe storms and operate safely over a period of several decades requires an assessment of their structural reliability against extreme wave loads and fatigue. Ocean waves and the responses of marine structures are usually modelled as stationary Gaussian random processes represented by stochastic Fourier series based on the work of Rice (1944, 1945). The standard stochastic model of a seastate assumes that the free-surface elevation may be represented by the superposition of a large number of Fourier components, typically 20–50 depending on the application, with random and independent phases drawn from the uniform distribution and amplitudes determined from the power spectral density of the seastate (Ochi 1998). The use of a seastate representation with a large number of independent sources of uncertainty in nonlinear wave body interactions leads to a computational task that may become prohibitively expensive when the statistics of extreme loads and responses are necessary. Therefore, a representation of the seastate with the smallest possible number of independent sources of uncertainty would be very beneficial. When coupled with the theoretical models of the nonlinear wave loads, explicit expressions for the nonlinear probability density functions and threshold up-crossing rates become possible. The significance of a representation of linear and nonlinear random signals in terms of a small number of independent canonical components was underscored by Kac & Siegert (1947) and Wiener (1958).

Sines and cosines are but one of all possible orthogonal basis function sets that may be used for the representation of deterministic or stochastic signals in the form of a series with stochastic coefficients. Given the power spectral density of the signal, an optimal orthogonal set of basis functions exists that fits the signal with the minimum number of independent sources of uncertainty. This basis set follows from the spectral decomposition theorems of Loeve (1945) and Karhunen (1947) (see Basilevsky 1994). An eigenvalue–eigenfunction problem is formulated and cast in the form of an integral equation of the first kind with a kernel equal to the autocorrelation function of the signal. The solution of this integral equation allows the representation of the signal as the linear superposition of the eigenfunctions multiplied by independent random variables with variance being equal to the corresponding eigenvalues. A rapid decay of the eigenvalues indicates that a small number of independent components are sufficient for the stochastic representation of the signal. This is indeed the case for narrow-banded seastates for which less than three to five terms in the Karhunen–Loeve (KL) representation may provide sufficient accuracy.

The solution of the KL eigenvalue problem is a complex task, and analytical solutions are available only in a few special cases. Such a case was considered by Slepian & Pollack (1961) and involves the *sinc* function as the kernel in the KL integral equation. The eigenfunctions turn out to be the prolate spheroidal wave functions (PSWFs), which form an orthonormal basis set with properties analogous to those of Fourier series, Chebyschev and Legendre polynomials. A distinct advantage of PSWFs is that they can be maximally concentrated in a given interval via the selection of a controllable bandwidth parameter known as the Slepian frequency. Moreover, they have the attractive property that they are invariant to a finite and infinite Fourier transform. Their application in communication theory is discussed by Moore & Cada (2004) and their mathematical properties and computation are presented in Xiao *et al*. (2001).

Drawing upon the properties of PSWFs, a general solution of the KL integral equation is developed in the present article. The kernel is given by the autocorrelation function of a stationary seastate with a known spectrum. The KL integral equation of the first kind is reduced to a spectral matrix eigenvalue problem that accepts an efficient solution. For narrow banded seastates, a rapid decay of the eigenvalues is observed. The number of dominant components in the Karhunen–Loeve–Slepian (KLS) representation depends on the magnitude of a tunable bandwidth parameter known as the Slepian frequency. Its selection for specific applications is discussed and its magnitude is found to depend on the bandwidth of the spectrum and the time scale of the memory effects in the wave–body interaction problem under consideration.

The use of the KLS representation of a seastate is discussed in connection with a recent theoretical nonlinear wave load model developed by Sclavounos (2012), which permits the derivation of explicit expressions for the extreme statistics for use in design. The use of the KLS representation for the design of a critical wave episode in a seastate is also described. The extension of the present solution methodology to the joint-KL representation of correlated stochastic ocean wave elevations and wind speed profiles is finally addressed, and its use discussed for the analysis and design of fixed and floating offshore wind turbines and for the forecasting of ocean wave elevations, wind speeds and power output of wind farms.

## 2. Spectral decomposition of stochastic wave records

Consider a wave elevation record *ζ*(*X*,*Y*,*t*) from a stationary and ergodic seastate where (*X*,*Y*,*Z*) is an earth fixed coordinate system with *Z*=0 being the calm water surface and with the positive *z*-axis pointing upwards. The wave record measured at *X*=*Y* =0 has the same statistical properties as the wave record measured at other (*X*,*Y*) locations. For the sake of brevity, the following notation is adopted for the signal *ζ*(*t*)≡*ζ*(*t*,0,0). The (*X*,*Y*,*Z*) coordinates are introduced again later in the article when the space–time dependence of the velocity potential corresponding to the seastate is derived.

The autocorrelation function of the signal is given by
2.1The expectation in (2.1) is taken over an ensemble of records at a fixed time *t* or over time for a single stationary record. For a stationary ergodic seastate, the two definitions lead to an identical result. In practice, the use of the temporal average is more convenient. It is hereafter assumed that *R*(*τ*) is a known deterministic even function of the time lag *τ*. The two-sided power spectral density of the signal *ζ*(*t*) is the Fourier transform of the autocorrelation function. The following standard relations hold:
2.2In (2.2), the standard deviation of the signal *σ* is also defined in terms of the one-sided wave spectrum *S*_{ζ}(*ω*), which accepts a number of standard parametrizations in ocean engineering discussed in Ochi (1998).

Consider a signal over a finite time interval (−*T*,*T*). The KL theorem states that *ζ*(*t*) accepts the expansion
2.3Assuming that *ζ*(*t*) is a stochastic process, the coefficients *α*_{n} are independent random variables such that
2.4The deterministic functions *f*_{n}(*t*) are solutions of an eigenvalue problem cast in the form of an integral equation of the first kind with the autocorrelation function as its kernel,
2.5A number of observations follow from the KL stochastic series expansion (2.3)–(2.5).

— The independent random variables

*α*_{n}are Gaussian, if the signal*ζ*(*t*) is Gaussian, which is often the case with ocean waves. Yet the expansion applies to non-Gaussian signals.— The eigenfunctions

*f*_{n}(*t*) are even and odd functions for positive and negative values of their argument in the range (−*T*,*T*). They are obtained from the solution of the eigenvalue problem derived later.— The rate of decay of the eigenvalues

*κ*_{n}with increasing*n*suggests the number of terms that are sufficient to keep in the stochastic series expansion (2.3). If this number is small, the signal is governed by a small number of independent sources of uncertainty with statistical properties given by (2.3)–(2.5).— The basis set

*f*_{n}(*t*) is optimal in the sense that it allows the representation of the autocorrelation function with the minimum number of terms in the series in the third equation of (2.5).— The KL representation maximizes the Shannon entropy measure, which reveals the minimum number of terms that are sufficient for the representation of the variability of the signal.

These and other optimality properties of the KL representation are discussed in Basilevsky (1994).

The analytical determination of the eigenfunctions and eigenvalues corresponding to the autocorrelation function of the signal is a non-trivial task. A general solution that applies to ocean waves and other signals with known power spectral densities and autocorrelation functions is developed in §3.

## 3. Solution of the Karhunen–Loeve integral equation

The autocorrelation function of the signal introduced in (2.5) is seen to accept an expansion in terms of the eigenfunctions *f*_{n}(*t*), which are yet unknown. These eigenfunctions should ideally have a shape similar to that of the autocorrelation function. Figure 1 plots the autocorrelation function of two seastates with significant wave heights *H*_{s}≃4*σ* equal to 6 m and 10 m, respectively. The spectral density shown in Figure 2 is the two-parameter Bretschneider spectrum (Ochi 1998, p. 35), and the autocorrelation function is its Fourier transform defined by (2.2). The standard deviation *σ* of the seastate elevation is also defined in (2.2) as an integral of the spectrum.

Figure 3 plots the *sinc* function, which is the autocorrelation function of a band-limited signal with a power spectral density being constant for −*Ω*<*ω*<*Ω* and zero for |*ω*|>*Ω*. The *sinc* function shares strong similarities with the seastate autocorrelation function and will guide the derivation of the solution of (2.5) described later, David Slepian observed that the eigenvalue problem (2.5) accepts an explicit solution when its kernel is the *sinc* function. This follows from the identity
3.1The eigenfunctions are the PSWFs *ψ*_{n}(*t*,*c*) that satisfy the following relations (Xiao *et al*. 2001):
3.2The identity (3.1) holds for all arbitrary *ψ*_{n}(*t*,*c*) functions. The PSWFs *ψ*_{n}(*u*,*c*) form an orthonormal basis set in the interval (−1,1). They are symmetric functions of *u* for even values of *n* and anti-symmetric functions of *u* for odd values of *n*. The eigenvalues *λ*_{n}(*c*) are real for the even and purely imaginary for the odd PSWFs. The PSWFs depend parametrically on the bandwidth parameter *c*=*ΩT* known as the Slepian frequency, which controls the degree of their decay outside the interval (−1,1). The bandwidth parameter must be selected properly for each application—a topic discussed in more detail later. The proper selection of the bandwidth parameter enlarges the space of available PSWFs, introducing a flexibility not present in Fourier series.

The first equation in (3.2) shows that the PSWFs are invariant with respect to the finite Fourier transform. This Fourier invariance of the PSWFs proves very convenient in surface wave problems because it permits the transition from time- to frequency-domain representations using the same basis set. This property will be used later in the article to derive the velocity potential and flow kinematics of a seastate starting from the representation of the free-surface elevation as the stochastic signal. Another property of PSWFs is that they decay at infinity for a large value of their argument at a rate that depends on the value of the bandwidth parameter *c*. This permits their use as an ambient wave signal that starts and ends with a calm free surface. This is an attractive property in numerical simulations carried out over a time window (−*T*,*T*) of finite duration that circumvents the need for an artificial ramp up of the computation. The decay of PSWFs at infinity also permits their use with conservation theorems that make use of the assumption that the ambient seastate and the radiation–diffraction disturbances generated by floating bodies decay at infinity (Sclavounos 2012).

The PSWFs are used next to solve the KL integral equation for a general two-sided power spectral density. This solution is used in the next section to generate a representation of a seastate and its associated velocity potential in terms of a small number of independent sources of uncertainty.

The third equation in (3.2) is a key property that allows the representation of the complex exponential in terms of a series of products of PSWFs. Its introduction in the definition of the autocorrelation function allows the reduction of the kernel of the KL integral equation into a separable form, which opens the way for the derivation of a symmetric matrix eigenvalue problem,
3.3The third equation in (3.3) provides a series expansion of the autocorrelation function in terms of PSWFs and the function *ϕ*_{j}(*τ*). This series allows the conversion of the convolution kernel into a separable one that is tractable analytically. Introducing the third equation of (3.3) in the KL integral equation (2.5), we obtain
3.4Solving for the eigenfunctions *f*_{n}(*t*) on the right-hand side, we obtain
3.5The eigenfunctions are cast in the form of a series of PSWFs with yet unknown coefficients . Introducing this series expansion on both sides of (3.4), we obtain
3.6Both sides of (3.6) involve a series summation over the PSWFs *ψ*_{j}(*t*/*T*,*c*). Multiplying both sides with the PSWF *ψ*_{i}(*t*/*T*,*c*), integrating over the interval −*T*<*t*<*T* and invoking their orthogonality from (3.2), we obtain
3.7Invoking the definition of *ϕ*_{i}(*τ*) from (3.3), it follows that
3.8In (3.8), the asterisk denotes the complex conjugate of the eigenvalue when it is purely imaginary for odd PSWFs. Substituting (3.9) in (3.8), we obtain a symmetric matrix eigenvalue problem of the form
3.9The eigenvalues *κ*_{n} and corresponding eigenvectors are indexed by *n* and provide the solution of the KL integral equation. The eigenfunctions are determined from (3.5)

The eigenvalue problem (3.9) consists of even and odd parts that are independent. This follows upon inspection of the elements of the matrix *E*_{ki}. Its elements are non-zero only if (*i*,*k*) are either both even or both odd. If one is even and the other odd, the last integral in (3.8) is identically zero by virtue of the symmetry of the spectral density *S*(*ω*) with respect to its argument and the symmetry/anti-symmetry of the PSWFs of even/odd order, respectively. The even/odd eigenvalue problems correspond to the even/odd eigenfunctions *f*_{n}(*t*) by virtue of their definition in (3.5). The eigenvectors follow from the solution of (3.9) for the even/odd problems, respectively. The rate of decay of the eigenvalues *κ*_{n} is rapid and determines the number of terms that are sufficient to keep in the KL series representation of the signal.

The elements of the matrix *E*_{ki} decrease with increasing (*i*,*k*) because the magnitude of the last integral in (3.8) decreases as (*i*,*k*) increase owing to the oscillatory behaviour of the PSWFs as their order increases. Numerical tests were performed with the infinite matrix *E*_{ki} truncated with increasingly large values of its order until no further sensitivity was observed in the magnitude of the eigenvalues *κ*_{n} and their rate of decay towards zero, as reported later. The truncated eigenvalue problem (3.8) and (3.9) was found to be well conditioned, and its solution was carried out using a standard matrix package and a subroutine for the evaluation of the PSWFs and their eigenvalues. Numerical results are presented in §5.

## 4. Flow kinematics of the Karhunen–Loeve–Slepian seastate

The solution of the KLS integral equation leads to the representation of the free-surface elevation at *X*=*Y* =0 in the form
4.1The coefficients *α*_{n} are independent random variables with zero mean and variances *κ*_{n} given by the eigenvalues of the eigenvalue problem (3.9). The *n*th basis function *f*_{n}(*t*) is seen from (4.1) to be a linear superposition of PSWFs with coefficients being equal to the elements of the corresponding eigenvector. The series in the first equation of (4.1) includes both even and odd basis functions corresponding to even/odd PSWFs used in the second equation of (4.1) to represent the even and odd basis functions *f*_{n}(*t*), respectively.

The representation (4.1) of the free-surface elevation and the Fourier invariance of the PSWFs stated in (3.2) allows the direct representation of the velocity potential in the fluid domain corresponding to (4.1). Using (3.2), the second equation in (4.1) may be cast in the form
4.2The last equation in (4.2) expresses the free-surface elevation corresponding to the basis function *f*_{n}(*t*) in terms of a Fourier integral over the frequency range −*Ω*<*ω*<*Ω*. It is noted that the eigenvalues *λ*_{j} are real for even PSWFs and purely imaginary for the odd PSWFs; therefore, the right-hand side of the last equation in (4.2) is always real.

Assuming a unidirectional seastate propagating in the positive *X*-direction, and invoking the linear dispersion relation, the velocity potential corresponding to the free-surface elevation record *f*_{n}(*t*) takes the form
4.3The velocity potential *ϕ*_{n}(*X*,*Z*,*t*) is a linear superposition of plane progressive monochromatic wave components with amplitudes given by the second equation in (4.3). The last equation of (4.3) states that the velocity potential corresponding to each Fourier component is the complex conjugate of itself when the sign of the frequency is reversed. This property, combined with the fact that ‘the eigenvalues *λ*_{j} of even PSWFs are real, while the eigenvalues of odd PSWFs are purely imaginary’, ensures that the velocity potential given by the first equation in (4.3) is a real quantity.

The velocity potential of the KLS representation of a stochastic seastate with elevation given by (4.1) follows in the form 4.4This completes the derivation of the velocity potential of the seastate. The flow velocity follows by taking the gradient of (4.4) and the pressure follows from Bernoulli's equation.

## 5. Numerical results

The numerical solution of the eigenvalue problem (3.9) was carried out with Matlab. The elements of the matrix *E*_{ki} involve the evaluation of integrals of the seastate spectral density multiplied by the product of two spheroidal wave functions, and was carried out by quadrature. The integration is restricted to the range −*Ω*<*ω*<*Ω* with the cut-off frequency *Ω* selected to include the power spectral density band of pertinence to each application. The PSWFs and their eigenvalues *λ*_{j} were evaluated using the algorithms presented in Xiao *et al.* (2001).

A key attribute of the PSWFs as an orthonormal basis set is their parametric dependence on the parameter *c*=*ΩT*. The proper selection of *c* for each application leads to the selection of PSWFs that lead to a rapid decay of the eigenvalues of (3.9), the minimum number of terms in the series expansion (4.1) and the minimum number of sources of uncertainty. In the case of the two seastate spectra plotted in figure 2, the cut-off frequency was selected to be *Ω*≃0.6 rad s^{−1}.

The magnitude of the cut-off frequency *Ω* in problems involving wave–body interactions will depend on the rate of decay of the power spectral density of the floating body response under study. The selection of the time window *t*∈(−*T*,*T*) depends on the degree of persistence of free-surface memory effects affecting the response under study, a longer time window being necessary for responses with strong memory effects. A seastate represented by a KLS expansion would be used as input in a time-domain numerical simulation of the nonlinear loads or responses of ships, offshore structures and fixed or floating wind turbines. The time window of the numerical simulation is typically selected so that transients of the response have died out and do not affect the useful part of the force or response simulation record. The time scale of transients is of the same order as the time scale of the memory effects of the load or response under consideration. The memory effects of the diffraction wave disturbance generated by the legs of offshore structures with diameters of 10–20 m or the floaters of wind turbines with diameters of 5 m are weak. Consequently, it is often sufficient to select the time window (−*T*,*T*) to be as wide as two wave periods. For a seastate modal period of *T*≃12 *s*, *c*=*ΩT*≃0.6×12=7.2. The memory effects associated with the responses of ships on the other hand are stronger and persist longer since ships are large-volume structures. Therefore, for wave load and ship response simulations, the time window (−*T*,*T*) needs to be 5–10 wave periods wide depending on the type of ship response of interest. For the same seastate modal period and cut-off frequency, this leads to values of *c*=*ΩT*≃18–36. Consequently, the PSWFs appropriate for each application have different bandwidth parameters *c*.

Figures 4 and 5 plot the first four even and odd basis functions *f*_{n}(*t*) obtained from the solution of the eigenvalue problem (3.8) and (3.9). The power spectrum corresponds to a seastate with *H*_{s}=10 m and peak period *T*_{p}=13.6, and the value of the Slepian frequency is *c*=*ΩT*=10.

Figure 6 plots the autocorrelation function of the seastate and its expansion in terms of the even basis functions *f*_{n}(*t*) given by (2.5). It may be seen that two to three terms in the series (2.5) are sufficient for the accurate representation of the exact autocorrelation function. The rapid convergence of the KLS expansion is also shown in figure 7, which plots the relative error of the variance of the signal given in (2.5) in the form of a KLS series of squares of the basis functions. It is seen that a 6 per cent error is achieved with three terms in the series, and near-perfect accuracy is achieved with four terms. The rate of decay of the even and odd eigenvalues is illustrated in figure 8.

The corresponding results for a Slepian frequency *c*=20, the same cut-off frequency of *Ω*=0.6 rad s^{−1} and twice as large a time window are shown in figures 9–13. It is seen that the convergence of the KLS expansion is very fast and near-perfect convergence of the variance is achieved with six terms in the series. The rate of convergence of the eigenvalues for the case with *c*=50 corresponding to the same cut-off frequency and a time window of 100 s is illustrated in figure 14.

It may be observed from figures 8, 13 and 14 that consecutive even and odd eigenvalues are nearly equal but not identical. In the Fourier series representation of a stochastic signal, the random coefficients multiplying the even basis functions (cosines) and the odd basis functions (sines) have the same variance as a consequence of the Pythagorian identity that the sum of the squares of sines and cosines equals one. This identity does not hold for the even and odd basis functions in the KLS series, which are not periodic functions of time, a desirable property in simulations. Numerical tests using a Fourier series representation of the autocorrelation function reveal that up to 40–50 terms including sines and cosines may be needed for a comparably accurate approximation of the autocorrelation function, the seastate spectral density and its variance defined by equations (2.2). Moreover, Fourier series lack the flexibility introduced by the tuning parameter *c*, which may be used to select the optimal number of terms in the KLS representation for a particular application.

An alternative measure of the variability of the signal is the Shannon entropy measure or information criterion. In the case of the KLS expansion, it is given as the sum of the product of the eigenvalues times their natural logarithm. Figure 15 plots the information criterion as a function of the number of terms in the KLS series for three values of the Slepian frequency *c*. Convergence is achieved with *n*=3, 5 and 15 terms for *c*=10, 20 and 50, respectively. No more information about the variability of the signal is revealed by including more terms in the series.

## 6. Discussion

The KLS representation of a seastate derived in this article enables the representation of the ocean environment with the minimum number of independent sources of uncertainty by selecting a tunable bandwidth parameter *c* optimally for each application. This representation proves particularly useful for the representation of critical wave episodes in a seastate and for the study of the extreme statistics and threshold up-crossing rates of nonlinear loads and responses.

### (a) Critical wave episodes

In a number of ocean engineering applications, the derivation of the shape of a critical wave episode with a given amplitude is of interest. This is often referred to as the design wave, which may be used for the evaluation of the loads and responses of ships and offshore structures. Traditionally, this wave is assumed to be a sinusoid with a large and prescribed amplitude. Yet, this is often a conservative assumption, and a more realistic representation is necessary.

Consider an ambient seastate with *H*_{s}=10 m and *T*_{p}=13.6. The shape of the most probable design wave in such a seastate with a height *H*_{d}=20 m may be readily determined from the KLS expansion. Assume that this wave will be used for the design of an offshore structure over a time window corresponding to a Slepian frequency of *c*=20. Without loss of generality, the maximum wave height is assumed to occur at *t*=0. It follows from (4.1) that
6.1Four terms were kept in (6.1), leading to a 6 per cent error in the KLS representation of the seastate elevation. The equality (6.1) enforces a linear constraint on the independent random variables *α*_{n}, which are assumed to be Gaussian, with variances being equal to the eigenvalues *κ*_{n}. Admissible values of *α*_{n} therefore lie on a hyper plane in a four-dimensional space. The values of *α*_{n} that lead to the most probable design wave are the ones that maximize the joint Gaussian probability density function of, *α*_{n},*n*=0,…,3,
6.2In (6.2), *f*(*x*) denotes the Gaussian probability density function, and the independence of the *α*_{n} has been invoked. The numerical solution of (6.2) for the values *α*_{n} that maximize their joint density subject to the constraint (6.1) is a straightforward maximization problem. By substituting the values of *α*_{n} obtained from the solution of (6.2) in (4.1), the shape of the design wave as a function of time is obtained. By definition, it achieves its maximum value *H*_{d} at *t*=0 and decays to zero for *t*<0 and *t*>0 when the even basis functions *f*_{n}(*t*) are expressed in terms of PSWFs as in (4.1). In this particular example, the odd basis functions do not enter the definition of the design wave because their values vanish identically at *t*=0.

A similar approach may be used for other design studies where instead of the design wave, other critical load episodes of interest in offshore applications are desired. The maximum-likelihood estimation described earlier for the selection of the coefficients *α*_{n} forms the basis of the first-order and second-order reliability methods in structural design discussed in Kiureghian (2000) for linear and nonlinear problems.

### (b) Nonlinear statistics of loads and responses

In most applications of interest in offshore engineering, the ambient seastate is assumed to be Gaussian. This is however not always an accurate assumption for the loads and responses of ships and offshore structures in severe seastates that exhibit significant departures from Gaussianeity. Nonlinearities in the loads and responses affect their probability distributions, which in turn affect significantly their extremes, threshold up-crossing rates and fatigue loads.

The literature on theoretical and computational models for the prediction of the nonlinear loads and responses of floating structures is extensive. Fully nonlinear simulations often lead to prohibitive computational efforts when used in a stochastic seastate for the estimation of the statistics of extremes via Monte Carlo simulations. Therefore, approximate methods are often used. Offshore structures typically consist of slender cylindrical members with diameters being small compared with the wavelength of waves encountered in severe seastates. Ships are slender, floating structures with beam and draft being small compared with the wavelength of design waves. These geometric attributes justify the linearization of the radiation–diffraction wave disturbance due to the body about the ambient seastate profile, which may be assumed to be nonlinear and with amplitude comparable to the characteristic dimension of the structure. This approximation was used by Sclavounos (2012) to derive compact expressions for the nonlinear time-domain loads on ships and offshore structures that are explicit functions of the ambient seastate kinematics. Combined with a parsimonious KLS representation of the ambient seastate, the stochastic nonlinear loads and responses may be expressed as quadratic and cubic functionals of a small number of random coefficients *α*_{n} that are independent Gaussian random variables with variance *κ*_{n}. The probability density functions of quadratic or cubic functionals of the random variables *α*_{n} assumed Gaussian are available in closed form in terms of their characteristic functions discussed by Kac & Siegert (1947) and Ochi (1998). This allows the derivation of explicit expressions of threshold up-crossing rates of extremes and fatigue loads for use in structural design discussed by Winterstein (1988).

An alternative approach often used for the study of the statistics of extremes is based on the derivation of the envelope of the process under consideration. In the case of the KLS representation of the free-surface elevation given by (4.1), the envelope is obtained as the square root of the sum of the squares of the free-surface elevation, and its conjugate process is defined as the Hilbert transform of the free-surface elevation. The Hilbert transform of the basis functions *f*_{n}(*t*) is a linear combination of the Hilbert transform of the PSWFs. The invariance of PSWFs with respect to the Fourier transform and an application of the Fourier convolution theorem yields the Hilbert transform of PSWFs in the form of a finite Fourier integral, which may be evaluated easily by quadrature. The probability distribution of the envelope may again be evaluated explicitly using the theory of characteristic functions.

### (c) Joint extreme statistics of wave and wind loads on offshore wind turbines

Offshore wind turbines based on fixed foundations or supported by floaters are exposed both to wind and wave loads. Their analysis and design therefore presents an increased level of complexity since models of the wind loads on the turbine must be coupled with models for the wave loads on the supporting structure. This task has been undertaken for floating wind turbines by Sclavounos *et al*. (2010) and Tsouroukdissian *et al*. (2011).

The design of offshore wind turbines also requires models for the joint extreme statistics of the wave and wind loads exerted on the structure in severe weather environments. The wind speed profile at hub height and the associated seastate elevation are correlated, and this correlation must be modelled and taken into account in the design of the offshore wind turbine. The number of sources of uncertainty governing the wind and wave environments is greater than those governing each individually. Therefore, the development of a joint-KL representation of the stochastic wind and wave loads that takes into account their correlation and leads to a minimum number of independent wind–wave sources of uncertainty would be very valuable.

Stochastic wind speeds are described by well-understood spectra, and experimental data exist that contain information about the covariance structure of wind speeds and seastate elevation. These attributes may be used to develop a joint wind–wave KL representation of a particular environment extending the methods developed in the present article. The autocorrelation function of wind speeds and their averages may be estimated from the Fourier transform of the Von Karman or Kaimal spectra that are routinely used by the wind energy industry. The rate of decay of the resulting autocorrelation function with increasing time lag is exponential like and conforms to a multi-step auto-regressive process. The KL representation of wind speeds may be carried out by starting from the KL expansion of the autocorrelation of wind speeds averaged over selected intervals. Such KL expansions are available in closed form for autocorrelation functions with exponential-like decay. Their explicit orthogonal eigenfunctions may be used as the basis function set for the derivation of the KL representation of wind speeds with more general but similar autocorrelation functions, in a manner analogous to the use of the *sinc* function for ocean waves described earlier. An analogous process may be followed for the derivation of the KL representation of the loads on wind turbine rotors obtained using standard computer programs.

Following the derivation of the KL representations of the wave loads on the fixed or floating foundation and the wind loads on the rotor, their statistical correlation may be studied by using the canonical correlation analysis discussed in Anderson (2003). The canonical correlation analysis leads to a parsimonious statistical model of the joint wave and wind loads on offshore wind turbines and their joint extreme statistics for use in design.

### (d) Forecasting of seastate elevation records and vessel responses

The ability to forecast wave elevation records and the responses of ships, offshore structures, floating wind turbines and wave energy converters has a wide range of applications for the safe operation of marine structures and the enhancement of the energy yield of ocean energy systems. Wave elevation records can be forecasted reliably several periods into the future in seastates that are reasonably narrow banded using multi-step autoregressive methods (Fusco & Ringwood 2010). The parsimonious KLS representation of ocean wave elevation records developed in the present article allows the derivation of explicit expressions of the minimum number of autoregressive coefficients necessary to predict the wave elevation record over future time intervals of interest for each particular application. This is the result of the optimality of the KL expansion and the selection of the optimal Slepian frequency. The coefficients of the auto-regressive model follows from the explicit solution of the Yule–Walker equations discussed in Box *et al*. (2008) using the representation of the KLS autocorrelation function in terms of PSWFs discussed earlier. In the present forecasting problem, knowledge of the signal is assumed over a finite time interval *T*. When the time interval *T* is infinite, the solution of the forecasting problem has been given by Wiener (1949). In the case of non-stationary seastates, the forecasting method described earlier may be generalized to account for time-dependent coefficients using a state space representation and Kalman filtering, as described in Durbin & Koopman (2001).

Applications of autoregressive forecasting models of seastate elevations and vessel responses include the prediction of ship responses in severe seastates in order to prevent excessive motions, ship-to-ship or ship-to-offshore-platform loading operation, the enhancement of the energy yield of compliant tension-leg platform floating wind turbines via the feed-forward control of the wind turbine blade pitch angle and the development of optimal algorithms for the maximization of the yield of wave energy converters.

## Acknowledgements

I thank Yeunwoo Cho and Thomas Luypaert for their programming assistance. This research has been supported by Enel, b_Tec and Alstom. Enel and b_Tec are members of the MIT Energy Initiative. This financial support is greatly appreciated.

- Received February 1, 2012.
- Accepted April 12, 2012.

- This journal is © 2012 The Royal Society