7.6 Coastal Processes

As waves enter shallow water, they undergo profound transformations — refraction, shoaling, breaking, and dissipation — that drive longshore and cross-shore sediment transport, shape coastlines, and create the dynamic surf zone environment. Understanding these processes is fundamental to coastal engineering and management.

Wave Refraction: Snell's Law

Waves approaching a coast at an angle slow down in shallower water, causing the crests to bend (refract) toward the depth contours. This is governed by Snell's law for waves:

$$\frac{\sin\alpha_1}{c_1} = \frac{\sin\alpha_2}{c_2} = \text{const along a ray}$$

where α is the wave angle relative to the bottom contours and $c = \sqrt{gh}$ in shallow water

Derivation: Snell's Law for Waves

Step 1: Wavefronts must be continuous across depth contours. Consider a crest line crossing from depth $h_1$ to $h_2$, making angles $\alpha_1$ and $\alpha_2$ with the contour normal.

Step 2: The time for a crest to travel between two adjacent rays must be equal on both sides of the depth boundary. Along the contour, the crest segment has length $\ell$. The distance traveled perpendicular to the crest is $\ell \sin\alpha$, covered at speed $c$:

$$\Delta t = \frac{\ell \sin\alpha_1}{c_1} = \frac{\ell \sin\alpha_2}{c_2}$$

Step 3: Cancel the common length $\ell$ to obtain the refraction law:

$$\frac{\sin\alpha_1}{c_1} = \frac{\sin\alpha_2}{c_2}$$

Step 4: This is identical to Snell's law in optics, where $c = c_0/n$ and $n$ is the refractive index. Here the "refractive index" is $n = c_0/c = c_0/\sqrt{gh}$, increasing as depth decreases — waves always bend toward shallower water, just as light bends toward denser media.

Refraction concentrates wave energy on headlands (convergent rays) and spreads it in bays (divergent rays), which is why headlands erode and bays accumulate sediment. The refraction coefficient $K_r = \sqrt{b_0/b}$ relates the spacing of wave rays to wave height changes. The combined shoaling and refraction coefficient determines the nearshore wave height:

$$H = H_0 \cdot K_s \cdot K_r, \quad K_s = \sqrt{\frac{c_{g0}}{c_g}}$$

where $K_s$ is the shoaling coefficient and $c_{g0}$ is the deep-water group velocity

Diffraction

Waves bend around obstacles (breakwaters, headlands) into the geometric shadow zone. Governed by the Sommerfeld diffraction solution. Critical for harbor wave agitation analysis.

Reflection

Vertical seawalls reflect waves (reflection coefficient R ~ 0.9–1.0). Rubble mound breakwaters absorb energy (R ~ 0.3–0.5). Reflected waves create standing wave patterns and scour.

Wave Breaking

Waves break when the particle velocity at the crest exceeds the phase speed. The breaking criterion relates wave height to water depth:

$$\gamma_b = \frac{H_b}{h_b} \approx 0.78$$

McCowan (1891) criterion; typical range 0.6–1.2 depending on beach slope and wave steepness

Derivation: Breaking Criterion

Step 1: In deep water, Stokes (1880) showed the limiting steepness for a stable wave is $H/L = 1/7$ (i.e., $H_{\max} \approx 0.142L$), where the crest forms a 120° angle. Beyond this, the crest becomes unstable.

Step 2: In shallow water, the wavelength is much larger than the depth, so steepness is no longer the limiting factor — instead the water depth limits the wave amplitude. The wave crest cannot rise higher than the available water column.

Step 3: McCowan (1894) solved for the maximum height of a solitary wave propagating in water of depth $h$, obtaining:

$$\frac{H_b}{h_b} = \frac{\pi}{4}\tanh^{-1}\left(\frac{\pi}{4}\right) \approx 0.78$$

Physical basis: Breaking occurs when the horizontal particle velocity at the crest equals the phase speed $c = \sqrt{gh}$. At this point, water at the crest outruns the wave form and the wave overturns. The ratio $\gamma_b$ varies from ~0.6 on gentle slopes to ~1.2 on steep slopes due to nonlinear effects.

The Iribarren number (surf similarity parameter) classifies breaker type:

$$\xi_0 = \frac{\tan\beta}{\sqrt{H_0/L_0}}$$

ξ < 0.5: Spilling

Gentle slopes. Foam cascades continuously. Wide surf zone. Most energy dissipation.

0.5 < ξ < 3.3: Plunging

Moderate slopes. Classic barrel/tube. Intense but localized dissipation. Surfable waves.

ξ > 3.3: Surging/Collapsing

Steep slopes. Wave surges up beach without breaking. Minimal dissipation. Strong backwash.

Radiation Stress and Longshore Current

Radiation stress is the excess momentum flux due to wave motion. Its gradients in the surf zone drive mean currents and water level changes:

$$S_{xx} = E\left(\frac{c_g}{c}\left(1 + \cos^2\alpha\right) - \frac{1}{2}\right), \quad S_{xy} = \frac{E c_g}{c}\sin\alpha\cos\alpha$$

Derivation: Radiation Stress

Step 1: Radiation stress is the depth-integrated excess momentum flux due to waves. The total onshore momentum flux through a vertical plane is:

$$S_{xx} = \int_{-h}^{\eta}\left(p + \rho u^2\right)dz - \int_{-h}^{0}\rho g(h+z)\,dz$$

where the second integral subtracts the hydrostatic (no-wave) contribution.

Step 2: For linear waves, substitute the pressure $p = -\rho gz + \rho g\eta \frac{\cosh k(z+h)}{\cosh kh}$ and horizontal velocity $u = a\omega \frac{\cosh k(z+h)}{\sinh kh}\cos\alpha$, then average over a wave period:

$$S_{xx} = E\left(n(1 + \cos^2\alpha) - \frac{1}{2}\right), \quad n = \frac{c_g}{c}$$

Step 3: In deep water $n = 1/2$, so $S_{xx} = E(\frac{1}{2}\cos^2\alpha)$. In shallow water $n = 1$, giving $S_{xx} = E(\frac{3}{2}\cos^2\alpha + \frac{1}{2})$. The gradient $\partial S_{xx}/\partial x$ drives wave setup, and $\partial S_{xy}/\partial x$ drives the longshore current.

Physical meaning: Radiation stress represents the flux of wave momentum. When waves shoal and break, momentum flux decreases, and this lost momentum is transferred to the mean flow — producing setup (raised water level) and longshore currents.

The longshore current driven by the oblique wave radiation stress gradient is:

$$V_l = \frac{5\pi}{16} u_m \sin\alpha_b \cos\alpha_b$$

where $u_m$ is the maximum orbital velocity at breaking and $\alpha_b$ is the breaker angle

Wave Setup/Setdown

Onshore radiation stress gradient raises mean water level in the surf zone (setup, up to ~20% of H_b) and lowers it seaward of breaking (setdown). Setup contributes significantly to storm surge.

Rip Currents

Narrow, seaward-directed jets driven by alongshore variations in wave height and setup. Speeds up to 2 m/s. Formed at gaps in sandbars, between headlands, or at structures. Dangerous to swimmers.

Sediment Transport and Beach Morphodynamics

The Shields parameter determines whether bed sediment is mobilized by the flow:

$$\tau^* = \frac{\tau_b}{(\rho_s - \rho)gd}$$

where $\tau_b$ is bed shear stress, $\rho_s$ is sediment density (~2650 kg/m³), d is grain diameter

Sediment motion initiates when $\tau^* > \tau^*_{cr} \approx 0.05$ (Shields critical value). The CERC formula estimates longshore sediment transport rate:

$$Q_l = \frac{K}{16(\rho_s/\rho - 1)(1 - p)} \sqrt{\frac{g}{\gamma_b}} H_b^{5/2} \sin(2\alpha_b)$$

K ≈ 0.39 (SI units), p is porosity (~0.4), $\gamma_b = H_b/h_b$

Derivation: CERC Formula

Step 1: Start from the longshore component of wave energy flux at the breaker line. The energy flux per unit crest length is $(Ec_g)_b$, and its longshore component is:

$$P_l = (Ec_g)_b \sin\alpha_b \cos\alpha_b = \frac{1}{2}(Ec_g)_b \sin(2\alpha_b)$$

Step 2: Empirically, the immersed-weight longshore transport rate $I_l$ is proportional to $P_l$: $I_l = K \cdot P_l$. Converting from immersed weight to volume: $Q_l = I_l / [(\rho_s - \rho)g(1-p)]$.

Step 3: Substitute $E = \frac{1}{8}\rho g H_b^2$ and the shallow-water group velocity $c_g \approx \sqrt{gh_b} = \sqrt{gH_b/\gamma_b}$:

$$Q_l = \frac{K}{16(\rho_s/\rho - 1)(1-p)}\sqrt{\frac{g}{\gamma_b}}\,H_b^{5/2}\sin(2\alpha_b)$$

Step 4: The empirical constant $K \approx 0.39$ (SI) was calibrated from field data (Komar & Inman, 1970). The $H_b^{5/2}$ dependence means transport is very sensitive to wave height — doubling $H_b$ increases transport by a factor of ~5.7.

Longshore Transport

Oblique waves drive sediment along the coast. Rate maximized at α_b = 45°. Typical rates: 10&sup4;–10&sup6; m³/yr. Creates spits, barrier islands, and drives shoreline evolution.

Cross-shore Transport

Storm waves (large, steep) move sand offshore into bars. Calm swell moves sand onshore, rebuilding beaches. Seasonal beach profiles: summer berm, winter bar.

Bedload vs Suspended Load

Coarse sediment (>0.5 mm) moves as bedload (rolling, saltating). Fine sediment (<0.2 mm) is carried in suspension. The Rouse number $P = w_s/(\kappa u_*)$ determines the transport mode: P < 1.2 suspended, P > 2.5 bedload.

Sediment Budget

A coastal sediment budget tallies all sources (river input, cliff erosion, onshore feed) and sinks (longshore transport out, offshore loss, dredging) for a littoral cell. The balance determines whether the shoreline advances, retreats, or remains stable.

Python: Wave Refraction, Longshore Current & Sediment Transport

Python

!/usr/bin/env python3

script.py130 lines

#!/usr/bin/env python3
"""coastal_processes.py - Refraction, longshore current, CERC formula"""
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import brentq

g = 9.81

def dispersion_k(omega, h):
    """Solve dispersion relation for wavenumber k."""
    k0 = omega**2 / g
    if k0 * h > 50:
        return k0
    func = lambda k: omega**2 - g * k * np.tanh(k * h)
    return brentq(func, 1e-8, max(k0 * 20, 10.0 / h))

# --- 1. Wave ray tracing over bathymetry ---
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# Simple bathymetry: parallel contours
nx_bath, ny_bath = 200, 300
x_bath = np.linspace(0, 2000, nx_bath)  # cross-shore (m)
y_bath = np.linspace(0, 3000, ny_bath)  # longshore (m)
X_b, Y_b = np.meshgrid(x_bath, y_bath)
h_bath = 0.02 * x_bath  # linear slope, 0 at shore, 40m at 2km
h_bath = np.maximum(h_bath, 0.5)

T = 10.0  # wave period
omega = 2 * np.pi / T
alpha0 = np.radians(30)  # initial approach angle

# Trace rays
nrays = 8
ax1 = axes[0, 0]
ax1.contourf(X_b, Y_b, np.tile(h_bath, (ny_bath, 1)),
             levels=20, cmap='Blues_r', alpha=0.5)

for ray_idx in range(nrays):
    y_start = 300 + ray_idx * 300
    x_pos, y_pos = [2000.0], [float(y_start)]
    alpha = alpha0

for step in range(500):
        h_local = max(0.02 * x_pos[-1], 0.5)
        k = dispersion_k(omega, h_local)
        c = omega / k

# Deep water values for Snell's law
        c0 = g * T / (2 * np.pi)
        sin_alpha = np.sin(alpha0) * c / c0
        sin_alpha = np.clip(sin_alpha, -1, 1)
        alpha = np.arcsin(sin_alpha)

ds = 10.0  # step size
        dx = -ds * np.cos(alpha)
        dy = -ds * np.sin(alpha)

new_x = x_pos[-1] + dx
        new_y = y_pos[-1] + dy
        if new_x < 10:
            break
        x_pos.append(new_x)
        y_pos.append(new_y)

ax1.plot(x_pos, y_pos, 'r-', lw=1.5)

ax1.set_xlabel('Cross-shore (m)'); ax1.set_ylabel('Longshore (m)')
ax1.set_title('Wave Ray Tracing (Refraction)')
ax1.set_xlim(0, 2000)

# --- 2. Longshore current profile ---
ax2 = axes[0, 1]
h_surf = np.linspace(0.5, 5, 100)
H_break = 2.0  # breaking wave height
alpha_b = np.radians(10)  # breaker angle
gamma_b = 0.78

# Longuet-Higgins longshore current model
x_norm = h_surf / (H_break / gamma_b)  # normalized position
V_max = (5*np.pi/16) * np.sqrt(g * H_break / gamma_b) * \
        np.sin(alpha_b) * np.cos(alpha_b)
# Simple profile: linear increase to breakpoint, decay beyond
V_long = np.where(x_norm <= 1.0,
                  V_max * x_norm,
                  V_max * np.exp(-2*(x_norm - 1)))

ax2.plot(h_surf, V_long, 'c-', lw=2)
ax2.axvline(H_break/gamma_b, color='r', ls='--', label='Break point')
ax2.set_xlabel('Depth (m)'); ax2.set_ylabel('Longshore current (m/s)')
ax2.set_title(f'Longshore Current (Hb={H_break}m, alpha={np.degrees(alpha_b):.0f} deg)')
ax2.legend(); ax2.grid(True, alpha=0.3)

# --- 3. CERC longshore transport ---
ax3 = axes[1, 0]
alpha_range = np.linspace(0, 90, 100)
K_cerc = 0.39
rho_s, rho = 2650, 1025
p = 0.4
Hb_vals = [1.0, 2.0, 3.0]

for Hb in Hb_vals:
    Ql = K_cerc / (16*(rho_s/rho - 1)*(1-p)) * \
         np.sqrt(g/gamma_b) * Hb**2.5 * \
         np.sin(2*np.radians(alpha_range))
    ax3.plot(alpha_range, Ql * 3600 * 24 * 365,
             lw=2, label=f'Hb = {Hb} m')

ax3.set_xlabel('Breaker angle (degrees)')
ax3.set_ylabel('Longshore transport (m^3/yr)')
ax3.set_title('CERC Longshore Sediment Transport')
ax3.legend(); ax3.grid(True, alpha=0.3)

# --- 4. Shields parameter ---
ax4 = axes[1, 1]
d_range = np.logspace(-4, -2, 100)  # grain size 0.1mm to 10mm
tau_b_vals = [0.5, 1.0, 2.0, 5.0]  # bed shear stress (Pa)

for tau in tau_b_vals:
    theta_s = tau / ((rho_s - rho) * g * d_range)
    ax4.loglog(d_range*1000, theta_s, lw=2, label=f'tau = {tau} Pa')

ax4.axhline(0.05, color='r', ls='--', label='Critical Shields')
ax4.set_xlabel('Grain size (mm)')
ax4.set_ylabel('Shields parameter')
ax4.set_title('Shields Parameter vs Grain Size')
ax4.legend(fontsize=8); ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('coastal_processes.png', dpi=150)
plt.show()

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Beach Morphodynamics and Equilibrium Profiles

Beach profiles reflect the balance between constructive (onshore) and destructive (offshore) sediment transport. The Dean equilibrium beach profile provides a first-order description:

$$h(y) = A y^{2/3}$$

where y is the offshore distance and A depends on grain size (A ≈ 0.1–0.2 m¹/³ for sand)

The Wright-Short classification defines six beach states from dissipative (flat, wide surf zone) to reflective (steep, narrow surf zone), parameterized by the dimensionless fall velocity:

$$\Omega = \frac{H_b}{w_s T}$$

where $w_s$ is sediment fall velocity; Ω < 1 reflective, 1 < Ω < 6 intermediate, Ω > 6 dissipative

Dissipative

Wide, flat surf zone with multiple bars. Spilling breakers. Fine sand. Typical of storm-dominated coasts (US East Coast). Low beach gradients (1:100).

Intermediate

Rhythmic bar and beach states (crescentic, transverse, ridge-runnel). Rip currents common. Most Australian beaches fall in this category. Complex 3D morphology.

Reflective

Steep beach face, no bar. Surging/collapsing breakers. Coarse sand or gravel. Low wave energy. Beach cusps common. Steep gradients (1:10 to 1:5).

Coastal Engineering Structures

Coastal structures modify wave energy and sediment transport. The Hudson formula provides the minimum armor unit weight for rubble mound breakwaters:

$$W = \frac{\rho_a H^3}{K_D (\rho_a/\rho_w - 1)^3 \cot\alpha}$$

where W is armor unit weight, H is design wave height, $K_D$ is stability coefficient, α is structure slope

Groynes

Perpendicular structures that trap longshore sediment, building up the updrift beach but starving the downdrift coast. Spacing typically 1–3 times groyne length for effective retention.

Beach Nourishment

Addition of compatible sand to the beach from offshore borrow areas. "Soft" engineering approach. Typical volumes: 100–500 m³/m of shoreline. Renourishment needed every 3–10 years depending on wave climate.

Detached Breakwaters

Offshore structures that reduce wave energy in their lee, creating salients or tombolos. Spacing and distance from shore control the morphological response. Often used in combination with beach nourishment for enhanced longevity.

Nature-Based Solutions

Mangrove forests, coral reefs, seagrass beds, and salt marshes provide natural wave attenuation (30–60% energy reduction). They adapt to sea level rise through vertical accretion and are increasingly recognized as cost-effective, self-maintaining coastal protection that also provides ecosystem services.

Fortran: 1D Cross-Shore Wave Propagation (Mild Slope Equation)

Fortran

1D cross-shore wave propagation using mild slope equation

program.f9099 lines

program mild_slope_equation
  ! 1D cross-shore wave propagation using mild slope equation
  ! with shoaling, refraction, and energy dissipation by breaking
  implicit none
  integer, parameter :: dp = selected_real_kind(15)
  real(dp), parameter :: pi = 3.141592653589793_dp
  real(dp), parameter :: g = 9.81_dp, rho = 1025.0_dp
  integer, parameter :: nx = 500
  real(dp) :: dx, T, omega, H0, alpha0
  real(dp) :: x(nx), h(nx), k_wave(nx), c(nx), cg(nx), n_ratio(nx)
  real(dp) :: alpha(nx), Ks(nx), Kr(nx), H(nx), E_flux(nx)
  real(dp) :: gamma_b, Hb, hb, V_long
  real(dp) :: k_local, c0, cg0, k0
  integer :: i, iter

T = 10.0_dp              ! wave period (s)
  omega = 2.0_dp * pi / T
  H0 = 2.0_dp              ! deep water wave height (m)
  alpha0 = 20.0_dp * pi / 180.0_dp  ! deep water approach angle
  gamma_b = 0.78_dp         ! breaking criterion
  dx = 5.0_dp               ! grid spacing (m)

! Deep water properties
  k0 = omega**2 / g
  c0 = omega / k0
  cg0 = c0 / 2.0_dp

! Setup cross-shore profile (planar beach, slope 1:50)
  do i = 1, nx
    x(i) = real(i-1, dp) * dx
    h(i) = max(0.02_dp * (real(nx-i, dp)) * dx, 0.3_dp)
  end do

! Solve dispersion relation at each point
  do i = 1, nx
    k_local = omega**2 / g  ! initial guess
    do iter = 1, 50
      k_local = omega**2 / (g * tanh(k_local * h(i)))
    end do
    k_wave(i) = k_local
    c(i) = omega / k_local
    n_ratio(i) = 0.5_dp * (1.0_dp + 2.0_dp*k_local*h(i) / &
                 sinh(2.0_dp*k_local*h(i)))
    cg(i) = n_ratio(i) * c(i)
  end do

! Snell's law for refraction angle
  do i = 1, nx
    alpha(i) = asin(min(1.0_dp, sin(alpha0) * c(i) / c0))
  end do

! Shoaling coefficient
  do i = 1, nx
    Ks(i) = sqrt(cg0 / cg(i))
  end do

! Refraction coefficient
  do i = 1, nx
    Kr(i) = sqrt(cos(alpha0) / max(cos(alpha(i)), 0.01_dp))
  end do

! Wave height with breaking dissipation
  do i = 1, nx
    H(i) = H0 * Ks(i) * Kr(i)
    ! Apply breaking limit
    if (H(i) > gamma_b * h(i)) then
      H(i) = gamma_b * h(i)
    end if
  end do

! Find breakpoint
  Hb = 0.0_dp; hb = 0.0_dp
  do i = 1, nx
    if (H0*Ks(i)*Kr(i) >= gamma_b * h(i)) then
      Hb = gamma_b * h(i)
      hb = h(i)
      exit
    end if
  end do

! Longshore current at breakpoint
  V_long = (5.0_dp*pi/16.0_dp) * sqrt(g*Hb) * &
           sin(alpha(i)) * cos(alpha(i))

! Output results
  print *, '=== Cross-Shore Wave Propagation ==='
  print '(A,F6.1,A,F6.1,A)', ' T = ', T, ' s, H0 = ', H0, ' m'
  print '(A,F6.1,A)', ' alpha0 = ', alpha0*180/pi, ' deg'
  print '(A,F6.2,A,F6.1,A)', ' Hb = ', Hb, ' m at depth = ', hb, ' m'
  print '(A,F6.3,A)', ' Longshore current Vl = ', V_long, ' m/s'

print *, ''
  print *, '  x(m)     h(m)    H(m)    alpha(deg)   Ks      Kr'
  do i = 1, nx, nx/20
    print '(F8.1, F8.2, F8.3, F10.2, F8.3, F8.3)', &
          x(i), h(i), H(i), alpha(i)*180.0_dp/pi, Ks(i), Kr(i)
  end do

end program mild_slope_equation

Click Run to execute the Fortran code

Code will be compiled with gfortran and executed on the server

← 7.5 Tsunamis Part 8: Climate and Oceans →