Classification: Internal Status: Draft

z
ξ(x,y,t)

x
w&

w
u u&

h Particle
position

From deterministic regular waves to a random field

Sverre Haver,
StatoilHydro,
January 2008

From a determinstic regular wave to a deterministic irregular

solution
Solution of linearized governing equations (i.e. boundary equations are applied at the mean
free surface and terms of order higher than 1 are neglected):

ω
ξ ( x, t ) = φ 0 sin(ωt − kx) = ξ 0 sin(ωt − kx)
g

Real ocean waves – except swell – do not

look like this. They look less regular!

Let us assume that the sea surface repeat ∞

⎡ 2 nπ 2 nπ ⎤
after T and use Fourier. Then we can write ξ (t ) = ∑ ⎢a n cos t + bn sin t⎥
the observed surface as sum of sinusoidals. n =1 ⎣ T T ⎦

Each component is a solution to the linearized equations, i.e. the sum is also a solution to the
linearized equations. Î we have an solution that looks more reasonable.

1
 3

Note:
Frequency resolution
A random solution is defined by length
of time window, T.
⎛b ⎞
By introducing: ξ n = a n2 + bn2 and θ n = arctan⎜⎜ n ⎟⎟
⎝ an ⎠
∞
2π
The sum can be written: ξ (t ) = ∑ξ
n =1
n cos(ω n t − θ n ) , frequency resolution (ωi+1 – ωi) is: Δω =
T
This solution is a determinstic solution. How could we establish a more general solution?

Let us introduce the phase as a random variable being uniformly distributed between 0 and 2π,
while keeping our original values for amplitude ξn, i.e.:
∞
Ξ (t ) = ∑ξ
n =1
n cos(ω n t − Θ n )

In order to obtain one realization, of Ξ(t), we generate N different phases, θn, n=1, N. By generating
a new set of N phases we will obtain a new realization of the random process, Ξ(t) a.s.o..

ÎWe have a solution from which we can generate an infinite amount of possible realizations. The
observed time history we started with is one of the many possible realizations.

Probabilistic structure of our solution I

∞
Solution: Ξ (t ) = ∑ξ
n =1
n cos(ω n t − Θ n )

What is the distribution of Ξ(t) = Ξt?

Observation: Ξt is a sum of many independent random components. If non of the components

dominates over the other, we will according to the central limit theorem have:

1 ⎧⎪ 1 ⎛ ξ ⎞ 2 ⎫⎪ Î Ξt is Gaussian with zero mean and

f Ξ (ξ ; t ) = exp⎨ − ⎜⎜ ⎟⎟ ⎬ variance equal to the sum of the component
2π σ Ξ (t ) ⎪⎩ 2 ⎝ σ Ξ (t ) ⎠ ⎪⎭
variances.

ξ n2
The variance of a harmonic wave is given by: σ Ξ2 =
n
2
N
The variance of Ξt thus becomes: σ (t ) = σ =
2
Ξ
2
Ξ ∑σ
n =1
2
Ξn

Observation: As far as we do not change the component amplitudes during the time window,
the process variance will be constant.

2
 5

Description of a Gaussian surface process

The joint distribution of Ξ(ti) and Ξ(tj) is a two-dimensional Gaussian:

1 ⎧⎪ 1 ⎡⎛ ξ ⎞
2
⎛ξ ⎞⎛ ξ j ⎞ ⎛ ξj ⎞
2
⎤ ⎫⎪
f Ξ i Ξ j (ξ i , ξ j ) = exp ⎨− ⎢⎜⎜ i ⎟⎟ − 2 ρ ij ⎜⎜ i ⎟⎟ ⎜⎜ ⎟⎟ + ⎜⎜ ⎟⎟ ⎥⎬
2π 1 − ρ σ ⎪⎩ 2(1 − ρ ij ) ⎣⎢⎝ σ Ξ ⎝σΞ ⎠⎝σ Ξ ⎠ ⎝σΞ ⎦⎥ ⎪⎭
2
2
ij
2
Ξ ⎠ ⎠

Where:

Cov (Ξ i , Ξ j )
ρ ij =
σ Ξ2
=
σ Ξ2
1
[ ]
E Ξ(t i ) Ξ (t j ) =
1
σ Ξ2
rΞΞ (τ ) , τ = ti - tj rΞΞ(τ) independent of absolute time, if
process is stationary.

A Gaussian stochastic process is in principle described by a joint distribution of its values at a large
number of times, t1, t2, …., tN. The parameters of this joint model can be determined from rΞΞ(τ).

Î A stationary Gaussian process is in a statistical sense completely characterized by the auto-

correlation function rΞΞ(τ).

Question: How can the variance be determined from the auto-correlation function?

Wave spectrum
In practise it is more common to work with the Fourier transform of auto-correlation function which is
the wave spectrum. The one-sided wave spectrum reads:

∞ ∞
1
∫e rΞΞ (τ ) = rΞΞ (−τ ) = ∫ e iωτ s ΞΞ (ω ) dω; τ ≥ 0
− iωτ
s ΞΞ (ω ) = rΞΞ (τ ) dτ ; ω ≥ 0
π −∞ 0

We know from question on previous slide that: σ Ξ2 = rΞΞ (0)

∞

From expression for the auto-correlation above, we see that: σ Ξ2 = ∫ sΞΞ (ω )dω Important
result!!
0

If we have expanded the process of length T in N component, the variance of the nth component
reads:
2π
σ Ξ2 = sΞΞ (ω n ) Δω
n
where Δω =
T

3
 7

Direct estimation of spectrum from a wave history

We know from a previous slide that the variance of the nth component also is given by:

ξ n2
σ Ξ2 =
n
2
Setting the two expressions for the variance equal, and solve with respect to the wave spectrum,
gives:

ξ n2
s ΞΞ (ω n ) =
2 Δω

The Fourier amplitides, ξn, are typically calculated from a time history of ξ(t), using a Fast Fourier
Transform technique.

Simulating a realization of X(t) corresponding to a known

wave spectrum
∞
The random process is given by: Ξ (t ) = ∑ξ
n =1
n cos(ω n t − Θ n )

From previous slide we have: ξ n = 2 s ΞΞ (ω n ) Δω

Î N
Ξ (t ) = ∑
n =1
2 s ΞΞ (ω n )Δω cos(ω n t − Θ n )

If we now generate m sets of N random phases, we get m realizations of the process. They will all
look different, but wave wave heights and wave periods will of course show some similarity. If we
calculate the wave spectrum from these simulated time histories, we will find that they all produce
(as expected) the input spectrum, i.e. a smooth spectrum if input spectrum is smooth.

If we estimate the wave spectrum from an observed time history of length T, we will find that the
estimated spectrum is of very irregular shape.
If our simulations shall reflect the same behaviour, we will also have to introduce the Fourier
amplitude, Ξn, as a random variable. In that case it should be introduced as a Rayleigh variable.

4
 9

Standard models for wave spectrum

J
Fully developed wind sea: Pierson-Moskowitz (PM)

sΞΞ ( f ) = 0.3125 hs2 t p− 4 f −5

{
exp − 1.25 t p− 4 f −4
}
PM
Growing wind sea: JONSWAP spectrum (J)
⎧ ⎛ f − f p ⎞ ⎫⎪
2
⎪
⎛ f ⎞
−5
⎧ ⎛ f ⎞ ⎫⎪
−4 exp ⎨ − 0.5 ⎜ ⎟ ⎬
⎪ ⎜ f p σ ⎠⎟ ⎪
exp⎨ − 1.25 ⎜ ⎟ ⎬ (1 − 0.287 ln γ )γ ⎩
⎪
s ΞΞ ( f ) = 0.3125 hs2 t p ⎜ ⎟
⎝ ⎭
T
⎜f ⎟ ⎜f ⎟
⎝ p⎠ ⎪⎩ ⎝ p ⎠ ⎪⎭

Pure swell sea: JONSWAP spectrum with γ e.g. 10-15

Combined sea (Wind sea plus swell): Torsethaugen spectrum (T)

Torsethaugen spectrum = superposistion of two Jonswap type spectra

Wind Sea

If spectrum given in f and one would prefer 1 ω

s ΞΞ (ω ) = s ΞΞ ( f = ) Swell Sea
to have in w, one can find the latter by: 2π 2π

10

Torsethaugen versus JONSWAP

Hs = 10.0m, Tp =10.0s Hs = 10.0m, Tp =20s

300 250
Spectral Density (m**2*s)
Spectral Density (m**2*s)

250 Torsethaugen (2004) 200 Torsethaugen(2004)

JONSWAP JONSWAP
200
150
150
100
100
50
50
0
0
0 0.04 0.08 0.12 0.16 0.2
0 0.04 0.08 0.12 0.16 0.2
Frequency (Hz)
Frequency (Hz)

Hs = 10.0m, Tp =14.0s

250
Torsethaugen (2004)
Spectral Density (m**2*s)

200 JONSWAP

150

100

50

0
0 0.04 0.08 0.12 0.16 0.2
Frequency (Hz)

5
 11

Range of validity for spectral models

(Figure is meant to give a qualitative indication.)

100-year
One should know which contour for
sea states that are of hs and tp
major importance when
selecting spectral model.

12

Second order surface process

The Gaussian surface model is an approximation. It is in agreement with linear wave theory.

ξ (t ) = ξ G (t ) + ξ 2 (t )
A second order model is given by:

Comparison Gaussian process

and second order sum process. Gaussian process Second order correction
In principle there will also be a
⎡
[ ]⎤
difference frequency component. N N
ξ 2 (t ) = Re ⎢∑∑ c m c n hmn
+
e i (ω m + ωn ) t −
+ hmn e i (ω m − ω n ) t ⎥
⎣ m =1 n =1 ⎦

If problem under consideration is sensitive

to the detail surface elevation, a non-Gaussian
model should be adopted.

Surface elevation is straight forward, the

corresponding kinematics are more challenging.

6
 13

Crest heigts of Gaussian and Scond order processes

1
0.9
0.8

Cumulative probability
0.7
0.6
0.5
0.4 Rayleigh
0.3 Jahns and Wheeler
0.2 Forristall
0.1 Winterstein
0
12 16 20 24 28
Crest height (m )

Distribution function for the 3-hour maximum crest height for a sea state with
hs = 18.2m and tp = 17s (10**(-4) annual probability sea states.
(Forristall = Second order)