6.5 Ocean Modeling

Ocean general circulation models (OGCMs) solve the equations of fluid motion to simulate ocean currents, temperature, salinity, and biogeochemistry. They are essential tools for understanding climate.

Governing Equations

Primitive equations (Navier-Stokes + approximations):

\( \frac{Du}{Dt} - fv = -\frac{1}{\rho}\frac{\partial p}{\partial x} + A_H \nabla^2 u + A_V \frac{\partial^2 u}{\partial z^2} \)

Momentum (x-direction)

Approximations

Boussinesq, hydrostatic, incompressible, thin shell

Additional Equations

Continuity, tracer (T, S), equation of state

Model Types

OGCM (Ocean General Circulation Model)

Full 3D physics. MOM, NEMO, POP, ROMS. 0.1°-1° resolution.

Coupled Models (AOGCM)

Ocean + atmosphere + ice + land. CESM, GFDL, HadGEM. Climate projections.

Earth System Models (ESM)

Add biogeochemistry, carbon cycle. CESM2, GFDL-ESM. IPCC scenarios.

Model Resolution

~100 km

Climate models

0.25°

~25 km

Eddy-permitting

0.1°

~10 km

Eddy-resolving

Mesoscale eddies (50-200 km) require at least 0.25° resolution. Submesoscale (<10 km) requires 1 km.

Python: Shallow Water Model

#!/usr/bin/env python3
"""ocean_modeling.py - Simple shallow water model"""
import numpy as np
import matplotlib.pyplot as plt

def shallow_water_step(u, v, h, dx, dt, f=1e-4, g=9.81, H=1000):
    """
    One time step of linearized shallow water equations
    """
    # Gradients
    dh_dx = np.gradient(h, dx, axis=1)
    dh_dy = np.gradient(h, dx, axis=0)
    du_dx = np.gradient(u, dx, axis=1)
    dv_dy = np.gradient(v, dx, axis=0)

    # Update velocities (momentum equations)
    u_new = u + dt * (f * v - g * dh_dx)
    v_new = v + dt * (-f * u - g * dh_dy)

    # Update height (continuity)
    h_new = h - dt * H * (du_dx + dv_dy)

    return u_new, v_new, h_new

# Initialize
nx, ny = 50, 50
dx = 100e3  # 100 km
dt = 3600   # 1 hour

u = np.zeros((ny, nx))
v = np.zeros((ny, nx))
h = np.zeros((ny, nx))

# Initial perturbation (Gaussian)
x = np.arange(nx) * dx
y = np.arange(ny) * dx
X, Y = np.meshgrid(x, y)
h = 1.0 * np.exp(-((X - nx*dx/2)**2 + (Y - ny*dx/2)**2) / (500e3)**2)

# Run model
for step in range(100):
    u, v, h = shallow_water_step(u, v, h, dx, dt)

plt.figure(figsize=(10, 8))
plt.contourf(X/1e6, Y/1e6, h, levels=20, cmap='RdBu_r')
plt.colorbar(label='SSH anomaly (m)')
plt.xlabel('x (1000 km)'); plt.ylabel('y (1000 km)')
plt.title('Shallow Water Model: Rossby Wave Adjustment')