Free Surface Waves

Simcenter STAR-CCM+ provides the VOF Waves model to simulate surface gravity waves on the interface between a liquid and a gas, typically the atmosphere. Application examples for VOF waves are the investigation of flow around ships or wave impact on offshore structures.

The concept of a steadily-progressing periodic wave train is a convenient model that is used in coastal and ocean engineering applications to give fluid velocities, pressures, and surface elevations caused by waves. Assuming that the waves are propagating steadily without change (the steady wave problem), the wave train can be uniquely specified and solved in terms of three physical length scales: water depth, wave length, and wave height.

Flat Wave

A flat wave represents a calm plane of water.

First Order Wave

A first order wave is modeled with a first order approximation to the Stokes theory of waves. This approximation generates waves that have a regular periodic sinusoidal profile.

The equation for horizontal velocity is:

$${\upsilon }_h=a\omega \cos (\mathbf{\text{K}}\cdot \mathbf{\text{x}}-\omega t)e^{Kz}$$

(2642)

The equation for vertical velocity is:

$${\upsilon }_v=a\omega \sin (\mathbf{\text{K}}\cdot \mathbf{\text{x}}-\omega t)e^{Kz}$$

(2643)

The equation for surface elevation is:

$$\eta =a\cos (\mathbf{\text{K}}\cdot \mathbf{\text{x}}-\omega t)$$

(2644)

where:

• $a$ is the wave amplitude
• $\omega$ is the wave frequency
• $\mathbf{\text{K}}$ is the wave vector
• $K$ is the magnitude of the wave vector
• $z$ is the vertical distance from the mean water level.

The wave period $T$ is defined as:

$$T=\frac{2\pi }{\omega }$$

(2645)

The wavelength $\lambda$ is defined as:

$$\lambda =\frac{2\pi }{K}$$

(2646)

The dispersion relation (between wave period $T$ and wave length $\lambda$ ) for first order waves in finite water depth $d$ is:

$$T={\left[\frac{g}{2\pi \lambda }\tanh \left(\frac{2\pi d}{\lambda }\right)\right]}^{-1/2}$$

(2647)

Whereas for infinite water depth, the dispersion relation is:

$$\lambda =\frac{g{T}^2}{2\pi }$$

(2648)

The wave shape is independent of depth.

Fifth Order Wave

A fifth order wave is modeled with a fifth order approximation to the Stokes theory of waves. This wave more closely resembles a real wave than a wave that is generated by the first order method. The wave profile and the wave phase velocity depend on the water depth, wave height, and current.

The fifth order VOF waves are based on work by Fenton [125].

The Ursell number $U_R$ is defined as [122]:

$$U_R=\frac{H{\lambda }^2}{d^3}$$

(2649)

where $H$ , is the wave height $\lambda$ is the wavelength and $d$ is the depth of the water.

This wave theory is valid only for Ursell numbers less than 30.

Superposition Wave

A superposition wave is a linear superposition of different partial first order waves. It can be used to simulate more complex wave phenomena, such as a cross sea or spectral waves. A cross sea is a sea state with two wave systems traveling at oblique angles. This state can occur when water waves from one weather system continue despite a shift in wind.

Cnoidal Waves

A cnoidal wave is a nonlinear and exact periodic wave solution of the Korteweg–de Vries equation, which describes the propagation of waves over a flat bed. Korteweg & de Vries [129] obtained periodic solutions which they termed ”cnoidal” because the surface elevation is proportional to the square of the Jacobian elliptic function $\mathrm{c}\mathrm{n}\left(\right)$ . Cnoidal waves are used to describe surface gravity waves that have a long wavelength when compared to the fluid depth. The cnoidal solution displays the long flat troughs and narrow crests that are observed in natural shallow water waves. As the wavelength tends to infinity, the cnoidal solution describes a solitary wave.

The cnoidal solution presents results in terms of expansions in wave height. The parameter that is used is the wave height relative to the trough depth of fluid, $H/h$ , which is denoted by $ϵ$ .

The series is expressed as a power series in $ϵ/m$ rather than $ϵ$ , [124]. As $m$ can be less than 1, it is better to monitor the magnitude of $ϵ/m$ than to have a power series in $ϵ$ with coefficients which are polynomials in $1/m$ .

The full solution to third order is presented. This solution is more applicable to shorter and reduced-amplitude waves, where the parameter $m$ is less than 0.96.

The symbol $\mathrm{c}\mathrm{n}$ is used to denote $\mathrm{c}\mathrm{n}(\alpha X/h|m)=\mathrm{c}\mathrm{n}(\alpha (x-ct)/h\left|m\right)$ .

Surface elevation calculation

$$\frac{\eta }{h}=1+\left(\frac{ϵ}{m}\right)m{\mathrm{c}\mathrm{n}}^2+{\left(\frac{ϵ}{m}\right)}^2(-\frac{3}{4}{m}^2{\mathrm{c}\mathrm{n}}^2+\frac{3}{4}{m}^2{\mathrm{c}\mathrm{n}}^4)+{\left(\frac{ϵ}{m}\right)}^3\left(\right(-\frac{61}{80}{m}^2+\frac{111}{80}{m}^3){\mathrm{c}\mathrm{n}}^2+(\frac{61}{80}{m}^2-\frac{53}{20}{m}^3){\mathrm{c}\mathrm{n}}^4+\frac{101}{80}{m}^3{\mathrm{c}\mathrm{n}}^6)$$

(2650)

Horizontal fluid velocity in the frame of the wave calculation

Variation in the $x$ -direction is assumed to be relatively slow and can be expressed in terms of a scaled dimensionless variable $\alpha x/h$ . $\alpha$ is a small quantity which expresses the relative slowness of variation in the $x$ -direction, and $h$ is the trough depth.

The series are expressed in terms of ${\alpha }^2$ (specifically $\delta =4{\alpha }^2/3$ ) [126]. The results are accurate even for high waves.

$$\frac{U}{\sqrt{gh}}=-1+\delta (\frac{1}{2}-m+m{\mathrm{c}\mathrm{n}}^2)+{\delta }^2(-\frac{19}{40}+\frac{79}{40}m-\frac{79}{40}{m}^2+{\mathrm{c}\mathrm{n}}^2(-\frac{3}{2}m+3m^2)-m^2{\mathrm{c}\mathrm{n}}^4+{\left(\frac{Y}{h}\right)}^2(-\frac{3}{4}m+\frac{3}{4}{m}^2+{\mathrm{c}\mathrm{n}}^2(\frac{3}{2}m-3m^2)+\frac{9}{4}{m}^2{\mathrm{c}\mathrm{n}}^4))+{\delta }^3(\frac{55}{112}-\frac{3471}{1120}m+\frac{7113}{1120}{m}^2-\frac{2371}{560}{m}^3+{\mathrm{c}\mathrm{n}}^2(\frac{71}{40}m-\frac{339}{40}{m}^2+\frac{339}{40}{m}^3)+{\mathrm{c}\mathrm{n}}^4(\frac{27}{10}{m}^2-\frac{27}{5}{m}^3)+\frac{6}{5}{m}^4{\mathrm{c}\mathrm{n}}^6+{\left(\frac{Y}{h}\right)}^2(\frac{9}{8}m-\frac{27}{8}{m}^2+\frac{9}{4}{m}^3+{\mathrm{c}\mathrm{n}}^2(-\frac{9}{4}m+\frac{27}{2}{m}^2-\frac{27}{2}{m}^3)+{\mathrm{c}\mathrm{n}}^4(-\frac{75}{8}{m}^2+\frac{75}{4}{m}^3)-\frac{15}{2}{m}^3{\mathrm{c}\mathrm{n}}^6)+{\left(\frac{Y}{h}\right)}^4(-\frac{3}{16}m+\frac{9}{16}{m}^2-\frac{3}{8}{m}^3+{\mathrm{c}\mathrm{n}}^2(\frac{3}{8}m-\frac{51}{16}{m}^2+\frac{51}{16}{m}^3)+{\mathrm{c}\mathrm{n}}^4(\frac{45}{16}{m}^2-\frac{45}{8}{m}^3)+\frac{45}{16}{m}^3{\mathrm{c}\mathrm{n}}^6))$$

(2651)

The leading term $-1$ arises because the wave is considered to be traveling in the positive $x$ -direction. The fluid is flowing under the wave in the negative $x$ -direction, relative to the wave, with velocities of the order of the wave speed.

Vertical fluid velocity calculation

The vertical fluid velocity can be obtained from Eqn. (329) by using the mass conservation equation $\partial U/\partial X+\partial V/\partial Y=0$ , and $d\left(\mathrm{c}\mathrm{n}\left(\theta \right|m)\right)/d\theta =-\mathrm{s}\mathrm{n}\left(\theta \right|m)\mathrm{d}\mathrm{n}\left(\theta \right|m)$ . Each term that contains ${(Y/h)}^i{\mathrm{c}\mathrm{n}}^j(\alpha X/h|m)$ , for $j>0$ , is replaced by $\alpha \mathrm{s}\mathrm{n}\left(\right)\mathrm{d}\mathrm{n}\left(\right)\left(\frac{j}{i+1}\right)\times {(Y/h)}^{i+1}{\mathrm{c}\mathrm{n}}^{j-1}\left(\right)$ . Hence, if Eqn. (329) is written as:

$$\frac{U}{\sqrt{gh}}=-1+\sum _{i=1}^5{\delta }^i\sum _{j=0}^{i-1}{\left(\frac{Y}{h}\right)}^{2j}\sum _{k=0}^i{\mathrm{c}\mathrm{n}}^{2k}\left(\right){\mathrm{\Phi }}_{ijk}$$

(2652)

where each coefficient ${\mathrm{\Phi }}_{ijk}$ is a polynomial of degree $i$ in the parameter $m$ , then the vertical velocity component is:

$$\frac{V}{\sqrt{gh}}=2\alpha \mathrm{c}\mathrm{n}\left(\right)\mathrm{s}\mathrm{n}\left(\right)\mathrm{d}\mathrm{n}\left(\right)\sum _{i=1}^5{\delta }^i\sum _{j=0}^{i-1}{\left(\frac{Y}{h}\right)}^{2j+1}\sum _{k=1}^i{\mathrm{c}\mathrm{n}}^{2(k-1)}\left(\right)\frac{k}{2j+1}{\mathrm{\Phi }}_{ijk}$$

(2653)

Irregular Waves

An irregular wave models wind seas (a short-term sea state) with a wave spectrum, that is, the power spectral density function of the vertical sea surface displacement. An irregular wave can use either the Pierson-Moskowitz spectrum or the JONSWAP spectrum [123]. Both spectra describe wind sea conditions that often occur for the most severe sea states.

The Pierson-Moskowitz spectrum $S_{PM}\left(\omega \right)$ is defined as:

$$S_{PM}\left(\omega \right)=\frac{5}{16}\left(H_S^2{\omega }_p^4\right){\omega }^{-5}\text{exp}\left(-\frac{5}{4}{\left(\frac{\omega }{{\omega }_p}\right)}^{-4}\right)$$

(2654)

where ${\omega }_p=\left(2\pi \right)/\left(T_p\right)$ represents the angular spectral peak frequency.

The JONSWAP spectrum $S_J\left(\omega \right)$ is expressed as:

$$S_J\left(\omega \right)=A_{\gamma }{S_{PM}\left(\omega \right)\gamma }^{\text{exp}\left(-0.5{\left(\frac{\omega -{\omega }_p}{{\sigma \omega }_p}\right)}^2\right)}$$

(2655)

where:

• $S_{PM}\left(\omega \right)$ is the Pierson-Moskowitz spectrum
• $\gamma$ is the non-dimensional peak shape parameter
• $A_{\gamma }=1-0.287\text{ln}\left(\gamma \right)$ is a normalizing factor
• $\sigma$ is the spectral width parameter:

  ${\sigma }_a$ for $\omega \le {\omega }_p$

  ${\sigma }_b$ for $\omega >{\omega }_p$ .