In [ ]:
from IPython.display import HTML

Altura e direção de onda usando a lei de snell e batimetria realística


In [ ]:
import numpy as np
import numpy.ma as ma

# Constantes.
g = 9.81
twopi = 2 * np.pi
    
# Relação de dispersão (Resolvida pelo método de Newton-Raphson).
def dispersion_relationship(T, h, Ho, thetao):
    omega = twopi / T
    Lo = (g * T ** 2) / twopi
    # Começa pela relação de dispersão de águas profundas.
    k = omega / np.sqrt(g)
    # Vamos aproximando por incrementos de "f".
    f = g * k * np.tanh(k * h) - omega ** 2
    while np.abs(f.max()) > 1e-7:
        dfdk = (g * k * h * (1 / (np.cosh(k * h))) ** 2 + g * np.tanh(k * h))
        k = k - f / dfdk
        f = g * k * np.tanh(k * h) - omega ** 2

    # Com o número de onda resolvido podemos calcular:
    L = twopi / k
    C = omega / k
    Co = Lo / T
    G = 2 * k * h / np.sinh(2 * k * h)
    
    # Aqui é basicamente a lei de Snell.
    theta = np.rad2deg(np.arcsin(C / Co * np.sin(np.deg2rad(thetao))))
    Kr = np.sqrt(np.cos(np.deg2rad(thetao)) / np.cos(np.deg2rad(theta)))
    Ks = np.sqrt(1 / (1 + G) / np.tanh(k * h))
    H = Ho * Ks * Kr
    # Retorna Altura de onda e ângulo.
    return map(ma.asanyarray, (theta, H))

Problema: Determine o período da onda que começaria a "sentir" o fundo em 10 metros de coluna d'água.

O fundo começa a afetar quando a coluna d'água é $\frac{L}{2}$, logo, para qual período de onda o comprimento de onda é igual a 20 metros?

Note que a altura de onda, na teoria linear, não afeta os cálculos.

Maui Jaws


In [ ]:
import os
import iris

longitude = lambda cell: -156.2765 <= cell <= -156.256
latitude = lambda cell: 20.934 <= cell <= 20.9405
constraint = iris.Constraint(coord_values=dict(longitude=longitude,
                                               latitude=latitude))

fname = '../../data/himbsyn.bathy_reduce.v19.nc'

if not os.path.isfile(fname):
    bathy = iris.load_cube('../../data/himbsyn.bathy.v19.grd', constraint)
    iris.save(bathy, fname)
else:
    bathy = iris.load_cube(fname)
print(bathy)

In [ ]:
import cartopy.crs as ccrs
import iris.quickplot as qplt
import matplotlib.pyplot as plt
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER

def plot_bathy(data, cmap=plt.cm.gist_earth_r):
    fig, ax = plt.subplots(figsize=(8, 8), subplot_kw=dict(projection=ccrs.PlateCarree()))
    qplt.pcolormesh(data, cmap=cmap)
    gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,
                      linewidth=1.5, color='gray', alpha=0.5, linestyle='--')
    gl.xlabels_top = False
    gl.ylabels_right = False
    gl.xformatter = LONGITUDE_FORMATTER
    gl.yformatter = LATITUDE_FORMATTER
    return fig, ax

In [ ]:
# Model requires positive values for bathymetry.
h = -bathy.data
land = bathy.data.mask
metadata = bathy.metadata

In [ ]:
# Start at the deepest bathymetry and goes from there up to the coast.
h = np.flipud(h)
rows, cols = h.shape

theta, H = [], []
for row in range(rows):
    if row == 0:  # First.
        thetao, Ho = dispersion_relationship(T=10, h=h[row, :],
                                             Ho=3, thetao=45)
    else:
        thetao, Ho = dispersion_relationship(T=10, h=h[row, :],
                                             Ho=Ho, thetao=thetao)
    H.append(Ho)
    theta.append(thetao)
    
# Flip then back.
H = np.flipud(H)
theta = np.flipud(theta)

# Apply land mask.
H = ma.masked_array(H, land)
theta = ma.masked_array(theta, land)

# Make them iris cubes for quickplotting.
H = iris.cube.Cube(H)
H.metadata = bathy.metadata

theta = iris.cube.Cube(theta)
theta.metadata = bathy.metadata

In [ ]:
fig, ax = plot_bathy(bathy)
fig, ax = plot_bathy(H, cmap=plt.cm.YlGnBu)
ax.set_title('Wave height [m]')
fig, ax = plot_bathy(theta, cmap=plt.cm.YlGnBu_r)
ax.set_title('Angle [degrees]')

In [ ]:
# http://www.coastal.udel.edu/faculty/rad/wavetheory.html
HTML('<iframe src=http://www.coastal.udel.edu/ngs/ width=700 height=350></iframe>')