Second part of the chapter on Coriolis from Introduction to Geophysical Fluid Dynamics, 2nd Edition
Here I implemented mfiles parabolic.m (script that calls the animation) and coriolisanim.m (a function for the animation)
The animation shows the trajectory of a particle on a parapoloid in absolute axes and the imprint of the trajectory on the rotating platform.
In [2]:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import animation
from JSAnimation import IPython_display
Other Examples (Figure 2.8):
- $\alpha$, $\beta$, X0 = 0, 0, 5 $\rightarrow$ part a.
- $\alpha$, $\beta$, X0 = 1, 0, 5 $\rightarrow$ part b.
- $\alpha$, $\beta$, X0 = -1, 0, 5 $\rightarrow$ part c.
In [3]:
n = 101 # Number of time steps to show a full period.
frames = n // 1 + 1 # Last time step (to show only half (quarter) a period divide by 2 (4)).
alpha = 1.2 # Scaled azimuthal initial velocity V_0/(X_0 Omega).
beta = 0.9 # Scaled radial initial velocity U_0/(X_0 Omega).
X0 = 3 # Initial distance (less then 10 for plotting reasons)
In [4]:
def coriolisanim(n, frames, alpha, beta, X0):
"""Function that shows the trajectory of a particle
on a parapoloid in absolute axes and the imprint of
the trajectory on the rotating platform."""
twopi = 2 * np.pi
pio2 = np.pi / 2
# Non-animate.
k = np.arange(1, 102)
xc = 10 * np.cos(k * twopi / 100)
yc = 10 * np.sin(k * twopi / 100)
XAA, YAA, xt, yt = [], [], [], []
fig = plt.figure(figsize=(6, 6))
ax = plt.axes(xlim=(-12, 12), ylim=(-12, 12))
ax.plot([10.5, 12], [0, 0], color='grey', linestyle='-')
ax.plot([0, 0], [10.5, 12], color='grey')
ax.plot(xc, yc, color='black')
ax.set_axis_off()
# Animated.
l1, = ax.plot([], [], color='grey', linestyle='-')
p1, = ax.plot([], [], linestyle='none', marker='o', mfc='black', mec='white')
l2, = ax.plot([], [], color='grey', linestyle='-')
l3, = ax.plot([], [], linestyle='none', marker='.', mfc='green', mec='white')
l4, = ax.plot([], [], linestyle='none', marker='o', mfc='green', mec='red')
l5, = ax.plot([], [], color='red', linestyle='-')
l1.set_data([], [])
l2.set_data([], [])
l3.set_data([], [])
l4.set_data([], [])
l5.set_data([], [])
p1.set_data([], [])
def init():
return l1, l2, l3, l4, l5, p1
def animate(j):
ot = j * twopi / n # omega t
dot = twopi / n # omega delta t
# Four line segments to see rotation of plane (two segments fixed, two
# segments attached to rotating table.
x1 = 7 * np.cos(j * twopi / n), 10 * np.cos(j * twopi / n)
z1 = 7 * np.sin(j * twopi / n), 10 * np.sin(j * twopi / n)
x2 = 10 * np.cos(j * twopi / n + pio2), 7 * np.cos(j * twopi / n + pio2)
z2 = 10 * np.sin(j * twopi / n + pio2), 7 * np.sin(j * twopi / n + pio2)
# XXA and YAA solution in absolute axes as a function of initial
# conditions.
XA = X0 * np.cos(ot) + beta * X0 * np.sin(ot)
YA = alpha * X0 * np.sin(ot)
# Storage of trajectory in absolute axes.
XAA.append(XA)
YAA.append(YA)
# Rotation of trajectory inprint on rotating axis.
for i in range(0, j):
xn = xt[i] * np.cos(dot) - yt[i] * np.sin(dot)
yn = xt[i] * np.sin(dot) + yt[i] * np.cos(dot)
xt[i] = xn
yt[i] = yn
xt.append(XA)
yt.append(YA)
# Adding of new point to the trajectory.
l1.set_data(x1, z1)
l2.set_data(x2, z2)
l3.set_data(xt, yt)
l4.set_data(XA, YA)
l5.set_data(XAA, YAA)
p1.set_data(x1[1], z1[1])
return l1, l2, l3, l4, l5, p1
return animation.FuncAnimation(fig, animate, init_func=init,
frames=frames, interval=100)
coriolisanim(n=n, frames=frames, alpha=alpha, beta=beta, X0=X0)
Out[4]:
In [5]:
HTML(html)
Out[5]:
This post was written as an IPython notebook. It is available for download or as a static html.
python4oceanographers by Filipe Fernandes is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.
Based on a work at https://ocefpaf.github.io/.