The Shape of Rogue Waves
Back in 2020, I got an unexpected email from Al Osborne, a physics Professor at the University of Torino and researcher at the Office of Naval Research in the US. I discovered that he is one the preeminent world experts on rogue waves, the 50-100 foot monsters that can arise even in moderate sea conditions and sink ships. Here's an excerpt from a BBC documentary on rogue waves. He turned out to be a fan of my work, years ago, on theta functions, as they produce soliton-type solutions of the non-linear Schrödinger equation which are a possible model for such waves. I was doubly fascinated because a) this was something that my student Emma Previato had worked out for her thesis (cf. her paper: Duke Math. J., 1985) and b) I have done a fair bit of ocean sailing and am most curious about such waves. And after struggling with the literature, it dawned on me that this also fits in with my work on the infinite dimensional manifold of simple closed plane curves and the idea of shape spaces. Let me explain.
I. Nonlinear gravity waves
Like almost all physics, one begins by simplifying the problem! Water is incompressible, ok, so their velocity vector field has no divergence. But their theory gets truly messy and complicated by their vorticity, the curl of that vector field. Well, don't forget that vorticity is preserved along streamlines in the absence of any external force. And when water truly settles down, as it does from time to time, even in mid-ocean (I have seen this and swam in deep ocean water as flat as a pancake), then its velocity vector field is zero! So mostly ocean water can be modeled by curl free divergence free vector fields. Sure, the wind is an external force and shelving bottoms create external forces near shores but in deep water and ignoring the topmost layers being blown around, it is irresistible to assume the curl is zero too. Aha, harmonic functions now make their appearance.
Let's do the math. First a domain: assume \( z \) is the vertical dimension and we wish equations for the time varying surface of an ocean \( \Omega^{(t)} \) of infinite depth. Denote the ocean's surface by \(\Gamma^{(t)}\) and its equation by \(z=\eta(x,y,t)\), (excluding breaking waves whose tops outrun the troughs). Let \( \vec v(x,y,z,t) \) be the velocity vector of the water. The motion of the surface is given by a normal vector field on \(\Gamma^{(t)}\) which must equal the normal component of \(\vec v\): $$\frac{\partial\Gamma^{(t)}}{\partial t}(P,t) = \vec v(P,t) \cdot \vec N_{\Gamma^{(t)}}(P), \text{ or } \frac{\partial\eta}{\partial t}(P,t)=\vec v(P,t) \cdot \left(-\frac{\partial\eta}{\partial x},-\frac{\partial \eta}{\partial y},1\right)$$ Next, there potential \( \phi(x,y,z,t) \) on \( \bigcup_t \Omega^{(t)} \), harmonic in \((x,y,z)\) such that \(\vec v = \nabla \phi\). Euler's equation becomes now the definition of the pressure: $$ \frac{\partial \phi}{\partial t} + \frac12\| \nabla \phi \|^2 = -p - gz$$ where we take the density of water to be 1, and g to be the force of gravity on earth's surface. However, on the surface, p must equal the atmospheric pressure, which we can absorb into the normalization of z, hence set p at the surface to zero. We assume that, at the bottom of the ocean, \(\phi\text{ and }\nabla\phi \rightarrow 0, z \rightarrow -\infty, p \rightarrow +\infty\). Finally, for every simply connected domain, one has the Poisson kernel \(\mathcal{P}_\Omega\) that computes every harmonic function on the domain from its boundary values. For flat seas, for instance, the domain is the lower half space and the kernel is \(-z/2\pi(x^2+y^2+z^2)^{3/2}\). Thus we complete the set of equations for the evolution of gravity waves using: $$\frac{\partial \phi}{\partial t} = -\mathcal{P}_{\Omega^{(t)}} \ast \left(gz+\tfrac12\|\nabla\phi\|^2\right)\big |_{\Gamma^{(t)}}$$
The majority of work on gravity waves deals with "wave trains", waves which are independent of one of the horizontal coordinates, e.g. y, leaving (x,z). In this case, \(\Omega^{(t)}\) can be taken as a plane domain and harmonic functions are the real parts of complex analytic functions of x+iz. Their real and imaginary parts are conjugate harmonic functions that determine each other by an integral transform generalizing the Hilbert transform. But very few people use these equations. Instead, they start with the ansatz: $$\eta(x,t) = \text{Re}\left(A(x,z,t).e^{i(kx-\omega t)}\right)$$ where A is a slowly varying "complex wave envelope". Then, by discarding judiciously terms thought to be small, one derives the result that A satisfies the non-linear Schrödinger equation with coefficients expressed in terms of k,ω. The beauty of this is that one has explicit solutions of the non-linear Schrödinger equation arising from theta functions on Jacobians of algebraic curves that appear to produce "rogue waves", (cf. Osborne's book Nonlinear Ocean Waves and the Inverse Scattering Transform, 2010). But wouldn't it be more fun to avoid the ansatz?
II. Shape spaces
Starting from completely different questions and motivations, Peter Michor and I had been studying, since the early 2000s, the infinite dimensional manifolds formed by the totality of a large variety of geometric structures. For example, if you fix an ambient manifold and look at all its submanifolds of some type, then the totality of such submanifolds is itself a manifold, albeit a pretty big one. Following algebro-geometric traditions, we called these the differentiable Chow manifolds. Riemann himself had noted the existence of such manifolds in his famous Habilitation lecture. There are many other examples but to fix ideas, the prime example, the one that has given rise to the most work, is this: take the ambient space to be simply the plane and consider in it all simple closed plane curves, making this a manifold in its own right. What continues to amaze me is the huge diversity of the geometric properties of this one space in the many natural metrics that it carries. A caution: I have on purpose not said how smooth or how jagged the curves are that define points in this space. Because of this, we don't have literally one space. It's exactly like the linear situation for function spaces: for each metric, there are distinct completions and these nest in each other in complex ways. OK, we have the same in the nonlinear realm: many instantiations, all being completions in different metrics of the core set of \(C^\infty\) curves. And there are finite-dimensional "approximations" like the space of non-intersecting n-gons. I'll give three examples of Riemannian metrics on this space.
Let's denote this core space by 𝒮 and its members by Γ with interior Ω. Then, as above, for all Γ, let \(T_\Gamma\text{ and }N_\Gamma\) be their tangent and normal bundles in the plane. A section of the normal bundle \(a:s \mapsto a(s).\vec N_\Gamma(s)\) represents a tangent vector to 𝒮 at the point representing Γ. A Riemannian metric on 𝒮 is then defined by a quadratic norm on every such section. The simplest possible one is just the \(L^2\) metric \( \|a\|^2 = \int_\Gamma a(s)^2ds.\) where s is arc-length. The resulting Riemannian manifold is a strange bird indeed: the infimum of path lengths between any two points of 𝒮 is zero! Geometrically, what's happening is that the sectional curvatures are all non-negative and, at any point unbounded so that conjugate points are dense on geodesics. Visually, the intermediate curves can grow spikes that shorten the above distance along any path as much as you want.
To get metrics that behave more normally, the standard way is to use Sobolev-type metrics. The best way to do this is by viewing 𝒮 as a quotient by the diffeomorphism group of the plane, Diff(R2) by the subgroup of diffeomorphisms that map the unit circle to itself. The Lie algebra of this group is the vector space of smooth vector fields on the plane and one can put Sobolev norms on them component-wise: $$ \|\vec v\|^2_{Sob-n} = \int_{R^2} ((I-\Delta)^n \vec v \cdot \vec v)dxdy.$$ If one extends this norm to be one-sided invariant and takes cosets on the same side, you get a quotient metric on 𝒮 for which the map from Diff to 𝒮 is a submersion: the tangent bundle "upstairs" splits into a vertical part tangent to the cosets and a horizontal part that is the pull back of the tangent bundle "downstairs". This is an isometry between the quotient metric on 𝒮 and the restriction of the one-sided invariant metric on Diff. All geodesics on 𝒮 for this metric lift to horizontal geodesics on Diff. A simple way to understand this definition is: $$ \| a \|^2_{Sob-n} = \inf\left\{ \|\vec v\|^2_{Sob-n}, \vec v \text{ on } R^2 \big| \vec v\cdot \vec N_\Gamma (s) = a(s) \right\}$$ In the land of pseudo-differential operators, there is such an Ln for which \( \|a\|^2_{Sob-n} = \int _\Gamma (L_n(a)a.ds\). Here n need not be an integer but, in all cases, Ln has degree 2n-1. So long as n > 1, these manifolds behave well, having geodesics and curvature etc., just like finite dimensional manifolds. Michael Miller's group at Johns Hopkins has used the 3D version of these metrics extensively to analyze medical scans.
What makes these Sobolev metrics really great is that, because arise from a one-sided invariant metric on a group, the lifted geodesics conserve their "momentum", it is transported by the diffeomorphisms in the lifted geodesic, leading to very simple geodesic equations. The cotangent space to 𝒮 at Γ can be thought of as the space of 1-forms ω to R2, but given only along Γ and that kill the tangent space of Γ. The inverse \(L_{Sob-n}^{-1}\) defines a norm here which has degree 1-2n, i.e. it's given by an integral kernel. Upstairs, the kernel is just convolution with a modified Bessel function, namely \(\|\vec x\|^{n-1}K_{n-1}(\|\vec x\|) \) times a constant. As Darryl Holm pointed out to me, if n > 1, this is a continuous function at 0 so the completion of the cotangent bundle contains δ functions. This means we can set the momentum to a sum of delta functions on Γ and get ODEs for the resulting geodesics which may be thought of as a kind of soliton. Note that the metric on the cotangent bundle is always weaker than that on the tangent bundle.
The final example is given by the Weil-Petersson metric in a suitable model of the universal Teichmüller space. One starts with the Riemann mapping from a) the inside of the unit disk to the inside of Γ, call this φint, and b) from the outside to the outside, called φext. The latter can be normalized by asking that infinity is mapped to infinity and that the derivative there is positive real. Then \(\psi = \phi_\text{ext}^{-1} \circ \phi_\text{int}\big|_{S^1}\) is a diffeomorphism of the circle called the welding map and is unique up to composition on the right by an a conformal self-map of the unit disk, i.e. a Möbius map. It can be shown that this map \(\mathcal{S} \mapsto \text{Diff}(S^1)\) creates an isomorphism between 𝒮 mod translations and scaling and the group of smooth diffeomorphisms of S1 modulo right multiplication by the three-dimensional Möbius subgroup of Diff. Once again we have a one-sided invariant metric on Diff(S1) via the formula: $$\| a(\partial/\partial\theta)\|^2 = \int_{S^1}H (a''' - a').a.d\theta = \sum_{n>0}(n^3-n)|\hat{a}(n)|^2$$ where prime is the θ derivative and H is the Hilbert transform for periodic functions. This defines a homogeneous norm of Sobolev degree 3/2 on 𝒮. The dual metric is given by a simple explicit continuous kernel. This is the famous Weil-Petersson metric. It turns out to be Kähler-Einstein metric with all negative sectional curvatures. The Einstein property says that its sectional curvatures must be small enough to make the Ricci trace finite, so in some sense, I think it is nearly flat. I think it's a gem of a space.
Essentially all the material in this section is available on my website, especially the notes from some Pisa lectures click here
III. Zakharov's Hamiltonian
Returning to the notation of the first part, the clue to linking these ideas on shape spaces to gravity waves is to consider the kinetic energy \(\int_\Omega\tfrac12\|\nabla \phi\|^2\) as a metric on 𝒮. OK, not exactly 𝒮 but now curves z = η(x) which are suitably tame at infinity (near the real axis) bounding a 2D slice of oceanic domains with infinite depth below them. Call this 𝒮Z. We assume the domain has fixed volume, meaning the mean of η is zero. A tangent vector to 𝒮Z at Γ is a normal vector field \(a(s)\vec N_\Gamma(s)\) to Γ such that \(\int_\Gamma a(s)ds = 0\). The Neumann boundary problem then defines a unique harmonic function in the interior with a as its normal derivative along the Γ and that also goes to 0 at -∞. If \(K_\text{Neu}\) is the corresponding Neumann kernel for the domain Ω, the metric is: $$\|a\|^2_Z = \iint_\Omega \tfrac12 \|\nabla \phi \|^2 dx dz, \quad \phi = K_\text{Neu}*a.$$ Note that because φ is harmonic, the integral can be rewritten: $$\iint_\Omega \|\nabla \phi \|^2 = \iint_\Omega \text{div}(\phi.\nabla\phi) = \int_\Gamma \phi.\frac{\partial \phi}{\partial n} = \int_\Gamma \phi.a$$ Thus we can interpret φ/2 as being the dual 1-form of the tangent vector a. In the simplest case η ≡ 0, \(K_\text{Neu}(s,x+iz) = \tfrac{1}{\pi} \log|s-(x+iz)| \) , hence \( \|a\|^2_Z = \tfrac{1}{2\pi}\iint \! a(s).a(t) \log|s-t| dsdt\). This is exactly the Sobolev H-½ norm because its Fourier transform is \( \int |\hat{a}(\xi)|^2 d\xi/\xi.\) So we are doing the opposite to what we did strengthening the L2 norm via derivatives. Here we have a weaker norm on the tangent bundle whose dual is stronger than it is.
On the other hand, to regularize the situation, we have potential energy as well as kinetic energy. This means the gravity wave equation is not a simple geodesic flow but a Hamiltonian flow where the potential is added to the norm squared. This is V.E. Zakharov's beautiful discovery in his 1968 paper Stability of Periodic Waves on the Surface of a Deep Fluid in the Zhurnal Prikladnoi Mekhaniki. The idea is to identify the cotangent space T *𝒮 with pairs (Γ,φ), φ harmonic on Ω and going to zero at -∞, taking Γ and φ as canonical dual variables. The Hamiltonian now is \(H(\Gamma, \phi) = \iint_\Omega \left( \tfrac12 \|\nabla \phi \|^2 + g.z\right)dxdz\) where the z term, after subtracting an infinite constant, should be interpreted as \(\int_\Gamma (g\eta(x)^2/2)dx\). One then checks that, if we write δΓ = a, then $$\delta H = \iint_\Omega \langle \nabla \phi, \nabla \delta \phi \rangle + \int_\Gamma \left( \tfrac12\|\nabla \phi\|^2 + g\eta\right) \delta\Gamma.ds$$ Rewriting the first term the way we did above for the metric, we find \( \iint_\Omega \langle \nabla \phi, \nabla \delta \phi \rangle = \int_\Gamma a.\delta\phi \), and we see that the Hamiltonian equations are the same as the equations for gravity waves: \(\frac{\delta H}{\delta \phi} = a = \frac{\partial \Gamma}{\partial t}\) and \(-\frac{\delta H}{\delta \Gamma} = -(\tfrac12\scriptsize{\|\nabla \phi\|^2 + g\eta})\big|_\Gamma = \frac{\partial \phi}{\partial t}\big|_\Gamma\).
Can we compute with such a system of equations? A key point is that the Hamiltonian is conformally invariant, hence one can shift everything to the unit disk using the time varying conformal map from the unit disk to Ω. This has been worked out by Dyachenko et al: see A.I.Dyachenko, E.A.Kuznetsov, M.D.Spector and V.E.Zakharov, Analytical Description of the Free Surface Dynamics of an Ideal Fluid, Physics Letters A, vol. 221, 1996 and V.E.Zakharov, A.I.Dyachenko and O.E.Vasilyev, New Method of Numerical Simulation of a Non-stationary Potential Flow of Incompressible Fluid with a Free Surface, European J. of Mechanics B - Fluids, vol.21, 2002. I would suggest however that an easy way to do numerical experiments is to replace the infinitely deep 2D ocean with the interior of a simple closed curve, close to the unit disk, as in the shape section above, while making gravity into a central force field based at the origin. Then Fourier series can be used and simulations without changing coordinates might be possible. Finding the rogue wave solutions by this route is a fascinating challenge and might even be of use to the study of genuine ocean rogue waves.