Solution (source code)

= Solution

The leading <Stokes drift> comes from evaluating the oscillatory <velocity> at the displaced particle position. <Taylor expansion> gives
$$
u(X,Z,t)-u(x,z,t)=(X-x)u_x+(Z-z)u_z+O(a^3k^2\omega).
$$
Here $u_x=-ak\omega e^{kz}\sin\vartheta$ and $u_z=ak\omega e^{kz}\cos\vartheta$. Inserting the initial-position displacements from the preceding solution yields
$$
\boxed{u(X,Z,t)-u(x,z,t)=a^2k\omega e^{2kz}[1-\cos(\omega t)]+O(a^3k^2\omega).}
$$
In particular the instantaneous difference is zero at $t=0$, as it must be. The unaveraged constant equality in the PDF is incompatible with its initial labels. Averaging over a period $\mathcal T=2\pi/\omega$ removes the oscillatory term:
$$
\boxed{u_s\equiv\left\langle u(X,Z,t)-u(x,z,t)\right\rangle=a^2k\omega e^{2kz}.}
$$
Equivalently, using mean parcel labels and the purely oscillatory displacements gives $a^2k\omega e^{2kz}(\sin^2\vartheta+\cos^2\vartheta)$ at this order. This recovers the intended constant result. \b[The period-mean drift is in the wave-propagation direction and decreases as $e^{2kz}$.] The corresponding second-order vertical difference is $a^2k\omega e^{2kz}\sin(\omega t)$ for initial labels and has zero mean. The resulting horizontal displacement per cycle is $u_s\mathcal T$; closed first-order circles do not imply zero second-order transport.