3.4 Local Winds

#!/usr/bin/env python3
"""
local_winds.py - Model sea breeze circulation

Run: python3 local_winds.py
Requires: pip install numpy matplotlib
"""
import numpy as np
import matplotlib.pyplot as plt

# Simple 2D sea breeze model
# Based on differential heating between land and water

def sea_breeze_circulation(t_hour, land_sea_contrast=10):
    """
    Simple model of sea breeze circulation strength
    t_hour: time of day (hours)
    land_sea_contrast: max temperature difference (K)
    """
    # Land-sea temperature difference varies with time
    # Max around 2-3 PM, reverses at night
    phase = (t_hour - 14) * 2 * np.pi / 24
    delta_T = land_sea_contrast * np.cos(phase)
    return delta_T

# Create time array
hours = np.linspace(0, 24, 100)
delta_T = sea_breeze_circulation(hours)

# Create a simplified 2D circulation diagram
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# Temperature difference over 24 hours
ax1 = axes[0, 0]
ax1.plot(hours, delta_T, 'r-', linewidth=2)
ax1.axhline(y=0, color='gray', linestyle='--')
ax1.fill_between(hours, 0, delta_T, where=(delta_T > 0), alpha=0.3, color='red', label='Land warmer')
ax1.fill_between(hours, 0, delta_T, where=(delta_T < 0), alpha=0.3, color='blue', label='Sea warmer')
ax1.set_xlabel('Hour of Day', fontsize=12)
ax1.set_ylabel('T(land) - T(sea) (K)', fontsize=12)
ax1.set_title('Land-Sea Temperature Difference', fontsize=14)
ax1.set_xlim(0, 24)
ax1.set_xticks(np.arange(0, 25, 3))
ax1.legend()
ax1.grid(True, alpha=0.3)

# Sea breeze schematic (afternoon)
ax2 = axes[0, 1]
x = np.linspace(-50, 50, 100)
z = np.linspace(0, 2000, 50)
X, Z = np.meshgrid(x, z)

# Simplified stream function for sea breeze
# Flow from sea to land at surface, return aloft
psi = np.sin(np.pi * X / 50) * np.exp(-Z / 800) * np.sign(50 - np.abs(X))

# Velocity components
u = np.gradient(psi, z[1]-z[0], axis=0)  # horizontal
w = -np.gradient(psi, x[1]-x[0], axis=1)  # vertical

ax2.streamplot(X, Z, u, w, color='blue', density=1.5)
ax2.axvline(x=0, color='brown', linewidth=3, label='Coastline')
ax2.fill_between(x[x >= 0], 0, -100, color='brown', alpha=0.3)
ax2.fill_between(x[x < 0], 0, -100, color='blue', alpha=0.3)
ax2.set_xlabel('Distance from coast (km)', fontsize=12)
ax2.set_ylabel('Height (m)', fontsize=12)
ax2.set_title('Sea Breeze Circulation (Afternoon)', fontsize=14)
ax2.set_ylim(0, 2000)
ax2.text(-25, 1500, 'SEA', fontsize=14, ha='center', color='blue')
ax2.text(25, 1500, 'LAND', fontsize=14, ha='center', color='brown')
ax2.annotate('', xy=(15, 200), xytext=(-15, 200),
             arrowprops=dict(arrowstyle='->', color='green', lw=3))
ax2.text(0, 300, 'Surface flow', ha='center', color='green', fontsize=10)

# Mountain valley wind schematic
ax3 = axes[1, 0]
# Create mountain profile
x_mtn = np.linspace(-50, 50, 200)
mountain = 2000 * np.exp(-x_mtn**2 / 200)

ax3.fill_between(x_mtn, 0, mountain, color='brown', alpha=0.5)
ax3.plot(x_mtn, mountain, 'k-', linewidth=2)

# Daytime upslope arrows
for x0 in [-30, -20, 20, 30]:
    ax3.annotate('', xy=(x0, mountain[np.argmin(np.abs(x_mtn-x0))]+200),
                 xytext=(x0, 200),
                 arrowprops=dict(arrowstyle='->', color='red', lw=2))

ax3.set_xlabel('Distance (km)', fontsize=12)
ax3.set_ylabel('Elevation (m)', fontsize=12)
ax3.set_title('Anabatic (Upslope) Wind - Daytime', fontsize=14)
ax3.set_ylim(0, 3000)
ax3.text(0, 2500, 'Warm air rises', ha='center', color='red', fontsize=10)

# Foehn effect
ax4 = axes[1, 1]
x_foehn = np.linspace(0, 100, 200)
mtn_profile = 3000 * np.exp(-(x_foehn-50)**2 / 200)

# Temperature profile
T_windward = 20 - 6 * mtn_profile / 1000  # moist adiabatic up
T_lee = T_windward.copy()
peak_idx = np.argmax(mtn_profile)
# Lee side: dry adiabatic descent (faster warming)
for i in range(peak_idx, len(x_foehn)):
    dz = (mtn_profile[i-1] - mtn_profile[i]) / 1000
    T_lee[i] = T_lee[i-1] + 9.8 * dz

ax4.fill_between(x_foehn, 0, mtn_profile, color='brown', alpha=0.3)
ax4.plot(x_foehn, mtn_profile, 'k-', linewidth=2)
ax4_twin = ax4.twinx()
ax4_twin.plot(x_foehn, T_windward, 'b--', linewidth=2, label='Windward (moist)')
ax4_twin.plot(x_foehn, T_lee, 'r-', linewidth=2, label='Lee (Foehn)')
ax4_twin.set_ylabel('Temperature (°C)', color='red', fontsize=12)
ax4.set_xlabel('Distance (km)', fontsize=12)
ax4.set_ylabel('Elevation (m)', fontsize=12)
ax4.set_title('Foehn Effect', fontsize=14)
ax4_twin.legend(loc='upper right')
ax4.annotate('', xy=(70, 1500), xytext=(30, 1500),
             arrowprops=dict(arrowstyle='->', color='purple', lw=3))
ax4.text(50, 1700, 'Prevailing wind', ha='center', color='purple', fontsize=10)

plt.tight_layout()
plt.savefig('local_winds.png', dpi=150, bbox_inches='tight')
plt.show()

print("Local Wind Types:")
print("  Sea Breeze: 5-10 m/s, penetrates 20-50 km inland")
print("  Land Breeze: 2-5 m/s, weaker than sea breeze")
print("  Katabatic: can exceed 50 m/s over ice sheets")
print("  Foehn/Chinook: can raise temperature 20°C in minutes")