Scattering of swell by currents

Abstract The refraction of surface gravity waves by currents leads to spatial modulations in the wave field and, in particular, in the significant wave height. We examine this phenomenon in the case of waves scattered by a localised current feature, assuming (i) the smallness of the ratio between current velocity and wave group speed, and (ii) a swell-like, highly directional wave spectrum. We apply matched asymptotics to the equation governing the conservation of wave action in the four-dimensional position–wavenumber space. The resulting explicit formulas show that the modulations in wave action and significant wave height past the localised current are controlled by the vorticity of the current integrated along the primary direction of the swell. We assess the asymptotic predictions against numerical simulations using WAVEWATCH III for a Gaussian vortex. We also consider vortex dipoles to demonstrate the possibility of ‘vortex cloaking’ whereby certain currents have (asymptotically) no impact on the significant wave height. We discuss the role of the ratio of the two small parameters characterising assumptions (i) and (ii) above, and show that caustics are significant only for unrealistically large values of this ratio, corresponding to unrealistically narrow directional spectra.


Introduction
Surface gravity waves (SGWs) play a key role in the exchanges of energy, momentum and gases between the ocean and the atmosphere (Villas Bôas & Pizzo 2021).SGWs are forced by the wind and modulated by ocean currents through transport and refraction.Over the past few decades, several studies have explored the effects of ocean currents on SGWs.Early theoretical work focusses on the formation of freak waves and identifies refraction as a possible mechanism for the generation of large amplitude waves (White & Fornberg 1998;Heller et al. 2008;Dysthe et al. 2008).
Recent studies examine how meso-and submesoscale ocean variability, such as fronts, filaments and vortices, induces a corresponding variability in wave amplitudes (Ardhuin et al. 2017;Romero et al. 2017Romero et al. , 2020;;Villas Boâs et al. 2020;Vrećica et al. 2022).These studies often characterise the wave amplitudes using the significant wave height H s , defined as four times the standard deviation of the surface displacement.They find that wave-current interactions at horizontal scales ranging from 10 to 200 km drive spatial gradients of H s at similar scales.This indicates that air-sea fluxes might have spatial variability on these relatively small spatial scales.
One common approach to studying wave-current interactions is the use of ray tracing, often in its simplest form in which the kinematics of SGWs is tracked by solving the ray equations and ray density is used as a proxy for wave amplitude (e.g., Kenyon 1971;Mapp et al. 1985;Quilfen & Chapron 2019).While this simple form of ray tracing is a valuable tool for understanding wave refraction, it does not provide an accurate quantification of changes in wave amplitude, in particular changes in H s .This quantification requires to solve the conservation equation for the density of wave action in the four-dimensional position-wavenumber phase space.This is challenging especially for the wave spectra of realistic sea states, distributed in both wavenumber and direction, instead of the pure plane waves that are often considered (see Heller et al. 2008, however).It is possible to solve the action equation numerically, albeit at great computational cost, either by discretising the phase space or by sampling its full four-dimensionality with a large ensemble of rays.
This paper proposes a complementary approach.It develops an asymptotic solution of the wave action equation, leading to explicit formulas for the changes in action and H s induced by localised currents.Motivated by their ubiquity in the ocean, we focus on swell, that is, SGWs characterised by a spectrum that is narrow banded in both frequency (equivalently, wavenumber) and direction.We exploit the smallness of two parameters reflecting the narrowness of the spectrum and the weakness of the current relative to the wave speed.We approximate the wave action equation to leading order and solve it in closed form by integration along its characteristics (the approximate ray equations) by inspection.The formulas we obtain show that the changes in action and H s depend on the currents through a 'deflection function' ∆ given by the integral of the vorticity along the primary direction of wave propagation.We apply these formulas to simple flows -vortices and dipoles -and compare their predictions with the results of full integrations of the action conservation equation by a numerical wave model.
We formulate the problem, relate action and H s , and introduce a model spectrum for swell in §2.We detail our scaling assumptions and carry out the (matched) asymptotics treatment of the wave action equation in §3.We compare asymptotic and numerical results for vortices and dipoles in §4.For vortices, we consider four different parameter combinations representative of ocean swell.We consider dipoles with axis along and perpendicular to the direction of the swell to demonstrate the possibility of a vanishing deflection function ∆, leading to asymptotically negligible changes in H s , a phenomenon we refer to as 'vortex cloaking'.In §5 we explore two limiting regimes of scattering: a linear regime, corresponding to weak currents and/or swell with relatively large angular spread, in which the changes in H s are linear in the current velocity, and a caustic regime corresponding to strong currents and/or small angular spread.The caustic regime, in which the changes in H s are large and concentrated along caustic curves, arises only for parameters values that are outside the range of typical ocean values.We conclude with a summary of our findings and discuss prospects for future work on the spatial variability of H s in §6.

Formulation
We study the scattering problem sketched in figure 1. Deep-water SGWs, with small initial directional spreading and a well defined peak frequency (swell) impinge on a spatially compact coherent flow, such as an axisymmetric vortex or a dipole.

Action conservation equation
In figure 1 we illustrate the scattering problem by tracing rays through an axisymmetric vortex.We go beyond ray tracing, however, by using asymptotic methods to obtain approximate analytic solutions of the conservation equation for the wave action density A(x, k, t) in the four-dimensional positionwavenumber space (Komen et al. 1996;Janssen 2004).The action conservation equation (2.1) relies on the WKB assumption of spatial scale separation between waves and currents.In (2.1) ω(x, k) is the absolute frequency of deep-water SGWs (2.2) We consider deep-water waves so that in (2.2) the intrinsic frequency is σ(k) = √ gk, with k = |k|.The current velocity is taken to be horizontal and independent of time and depth, U (x) = U (x, y) x + V (x, y) ŷ. (2. 3) The wave action equation (2.1) provides a phase-averaged description of the scattering problem made possible by the scale separation between waves and currents.This places our work in contrast to that of Coste et al. (1999), Coste & Lund (1999) and McIntyre (2019) who examined scattering without the simplification afforded by scale separation and discuss phase effects such as the Aharonov-Bohm effect.We also assume fixed currents and do not consider how these might be modified by the presence of waves (see e.g.Humbert et al. 2017;McIntyre 2019).

Action spectrum and significant wave height
Denoting the sea-surface vertical displacement by ζ(x, t), with root mean square ζ rms , and following Komen et al. (1996), we introduce a spectrum F(k, x, t) such that Later we use a polar coordinate system (k, θ) in k-space so that in (2.4) dk = k dkdθ.The kinetic and potential energy densities for deep-water SGWs are equipartitioned so that the energy spectrum is gF and the action spectrum, A(x, k, t) in (2.1), is A = gF/σ.The significant wave height, 4ζ rms (Komen et al. 1996), is therefore (2.5) Figure 1: The scattering problem: a localised flow, here shown as an axisymmetric vortex with radius rv , scatters waves incident from the left (x → −∞) with action spectrum A⋆(K, Θ).Rays bend significantly only in the scattering region in which there is non-zero vorticity i.e.where x = O(rv ).In this illustration rv is equivalent to ℓs .(a) The case δ ̸ = 0: directional spreading in the incident spectrum A⋆ is indicated schematically by two rays emanating from each source point.(b) The case δ = 0 (or much less than ε): the incident spectrum A⋆ is a plane wave with little or no directional spreading.
The incident swell is characterized by a spatially uniform spectrum F ⋆ (k) with constant significant wave height H s⋆ .The subscript ⋆ denotes quantities associated with the incident waves.Swell is characterized by a narrow spectrum in both wavenumber k (equivalently, frequency σ) and direction θ.The dominant wavenumber of the incident swell is k ⋆ with frequency σ ⋆ = √ gk ⋆ , and the dominant direction is taken without loss of generality as θ = 0. Thus, as illustrated in figure 1, the waves arrive from x = −∞ and impinge on an isolated flow feature, centred at (x, y) = (0, 0).As an example of incident spectrum we use a separable construction described in appendix A. In the narrow-band limit corresponding to swell, this spectrum simplifies to the Gaussian . (2.6) The two parameters δ k and δ θ capture the wavenumber and directional spreading (see Appendix A).The narrow-band limit assumes that δ k /k ⋆ ≪ 1 and δ θ ≪ 1.

The scattering problem
We consider an incident spectrum such as (2.6).To make its localisation in k and θ explicit we introduce the O(1) independent variables Focus on Fluids articles must not exceed this page length where δ ≪ 1 is a small dimensionless parameter.The incident action spectrum has the form A(x, y, k, θ) = A ⋆ (K, Θ) as x → −∞, (3.2) where the function A ⋆ (K, Θ) is localised where both K and Θ are O(1).The example spectrum (2.6) is of this form provided that δ k /k ⋆ and δ θ are both O(δ).This assumption of similarly small spectral widths in k and θ enforces the relevant distinguished limit for the scattering problem.
We assume that the currents are weak (e.g.Peregrine 1976;Villas Bôas & Young 2020).This means that the typical speed U of the currents is much less than the intrinsic group velocity of the incident swell c ⋆ : (3.4) Accordingly we rewrite the frequency (2.2) as We indulge in a slight abuse of notation here: we develop the approximation in dimensional variables, hence the dimensionless parameters ε and δ in expressions such as (3.1) and (3.5) should be interpreted as bookkeeping parameters to be set to one at the end.We examine the distinguished limit and use matched asymptotics to solve the action conservation equation (2.1).We emphasise that γ = O(1) is a formal assumption that enables us to capture the broadest range of relative size of ε and δ, including ε ≪ δ and δ ≪ 1 (see §5).
3.1.The scattering region: x = O(ℓ s ) The spatially compact flow has a typical horizontal length scale which we denote by ℓ s .We refer to the region where x = O(ℓ s ) as the 'scattering region'.The solution in this region has the form A(K, Θ, x, y) (3.7) and must limit to With A in (3.7) the transport term in (2.1) is approximated as In particular, transport by the current, εU • ∇ x A is negligible compared with transport by the intrinsic group velocity c ⋆ .With the approximations (3.10) the refraction term in (2.1) simplifies to Thus in the scattering region the leading-order approximation to (2.1) is One might solve (3.12) using its characteristics -the ray equations -or by inspection.By either method the solution to (3.12) that matches the incident action spectrum (3.2) as x → −∞ is found to be It is insightful to introduce the vorticity Z def = V x − U y and write (3.13) as (3.14)For reference, we rewrite this expression in terms of the original independent variables, setting the bookkeeping parameters ε, δ, and hence γ to 1 to obtain where we have introduced the dimensionless 'deflection' (3.17) According to (3.16) the effect of the flow on the dependence of A on K is reversible: after passage through the scattering region this dependence reverts to the incident form.In contrast, there is a net change in Θ, quantified by the deflection ∆(y).This can be related to classical scattering of particles by viewing y as the impact parameter of a wavepacket.The scattering cross section, defined as dy/dθ ∞ where θ ∞ is the angle of propagation of the wavepacket as x → ∞, is then −1/(ε∆ ′ (y)).
To physically interpret (3.16) and ∆(y), recall that if ε is small then ray curvature ≈ vorticity group velocity , (3.18) The approximation in (3.18) requires only ε ≪ 1 (e.g.Kenyon 1971;Landau & Lifshitz 2013;Dysthe 2001;Gallet & Young 2014).Passing from (3.18) to (3.19) requires the further approximation that k is close to k ⋆ so that the group velocity in the denominator of (3.18) can be approximated by the constant c ⋆ .On the left of (3.18) ray curvature is dθ/dℓ, where ℓ is arc-length along a ray.But within the compact scattering region we approximate ℓ with x.Thus the deflection ∆(y) in (3.17) is the integrated ray curvature, accumulated as rays pass through the scattering region in which x = O(ℓ s ) and vorticity Z(x, y) is non-zero.From (3.17) and (3.18) we conclude that the scattering region is best characterized as the region with O(1) vorticity, e.g. the vortex core in figure 1 (hence ℓ s = r v with r v a typical vortex radius).The region with palpably non-zero velocity is much larger.In figure 1 the rays are straight where x = O(r v /ε), despite the slow (∝ r −1 ) decay of the azimuthal vortex velocity.
3.3.The far field: x = O(ℓ s /δ) Far from the scattering region, where x ≫ ℓ s , we introduce the slow coordinate X def = δx.In the far-field the currents and hence the refraction term ∇ x ω • ∇ k A in (2.1) are negligible.The steady action conservation equation collapses to (3.20) i.e. propagation along straight rays.Retaining only the leading-order term gives By inspection the solution of (3.21) that matches the intermediate solution (3.22)This formula, which converts the incident spectrum into the far-field spectrum, is a key result of the paper.In terms of the original independent variables and with the bookeeping parameters set to 1 it takes the convenient form (3.23) 3.4.Significant wave height Significant wave height H s is the most commonly reported statistic of wave amplitudes, being routinely observed by satellite altimeters and wave buoys.We obtain an approximation for H s by performing the k and θ integrals in (2.5) using the approximations (3.15) and (3.23) for A(x, k).
The scattering region is simple.We can approximate σ and dk in (2.5) by σ ⋆ = σ(k ⋆ ) and k ⋆ dkdθ to find The second equality holds because, according to (3.15), A(x, k) is obtained from A ⋆ (x, k) by an x-dependent shift of the k and θ that does not affect the integral.Thus H s in the scattering region is unchanged from the incident value H s⋆ .This conclusion also follows directly from steady-state wave action conservation under the assumptions ε, δ ≪ 1: multiplying (3.12) by σ ⋆ k ⋆ and integrating over k and θ we find Hence H s (x) = H s⋆ throughout the scattering region.
In the far field, H s is obtained by substituting (3.23) into (2.5).The result is The k-integral can be evaluated in terms of the incident directional spectrum which, in the general case of a non-separable spectrum, is defined as (3.28) We summarize the results above with: and velocity where r 2 = x 2 + y 2 .The vortex radius r v can be taken as the scattering length scale ℓ s .The maximum azimuthal velocity is U m = 0.072 κ/r v at radius 1.585 r v .The deflection (3.17) resulting from this Gaussian vortex is The asymptotic solution in the scattering region is obtained from (3.15) as where erf is the error function.Eq. (4.4) can be combined with the far-field approximation (3.23) into a single, uniformly valid approximation, The significant wave height is approximated by (3.29) which can be written as the uniform expression where x + is equal to x for x > 0 and to 0 for x < 0 and (4.3) is used for ∆.
We now compare the matched asymptotic (MA hereafter) predictions (4.5)-(4.6)with numerical solutions of the wave action equation ( 2 The parameter s > 0 controls the directional spreading: for s ≫ 1, (4.7) reduces to the Gaussian in (2.6) with directional spreading δ θ = 2/s.The configuration of WW3 and spectrum parameters are detailed in Appendix B. The most important parameter is the peak frequency of the incident spectrum, taken fixed for all simulations as σ ⋆ = 0.61 rad s −1 .This corresponds to a period of 10.3 s, wavelength of 166 m and group speed c ⋆ = 8 m s −1 .Because the problem is linear in the action density, the values of ζ rms⋆ or equivalently H s⋆ are less important.
For definiteness we set H s⋆ = 1 m. Figure 2 compares the wavenumber-integrated wave action A(x, y, k, θ) dk obtained from (4.5) and WW3 for a Gaussian vortex with maximum velocity U m = 0.8 m s −1 and directional spreading parameter s = 40.Figure 2 shows a good agreement, especially in the far-field region (x ⩾ 3r v ).The most noticeable difference between MA and WW3 is in panels c and d, which show a section through the middle of the vortex.The MA action spectrum in panel d is obtained via a y-dependent shift in A ⋆ (k, θ); there is no change in the intensity of A associated with this shift.In panel c, on the other hand, the intensity of the WW3 action spectrum varies with y/r v .We attribute this difference to asymptotically small effects such as the contribution U • ∇ x A to wave-action transport.
In the remainder of this section, we assess the dependence of significant wave height H s on the directional spreading parameter s and flow strength U m .We consider the four different combinations of s and U m given in Table 1.The corresponding values of the dimensionless parameters, taken as are also in the table.
Observations of the directional spreading for swell typically range between 10 • − 20 • (Ewans 2002), which correspond to a range of s between 16 and 66.In our experiments, setting s = 10 and s = 40 leads to a directional spreading of 24 • and 12 • respectively, which correspond to very broad and very narrow swells.Figures 3 and 4 show the significant wave height anomaly for each combination of s and U m .Because of our choice of H s⋆ = 1m, h s in cm can be interpreted as the fractional change in significant wave height expressed as a percentage.A control run of WW3 in the absence of currents shows that h s is not exactly zero but decreases slowly with x.This is caused by the finite yextent of the computational domain which leads to a wave forcing with compact support.To mitigate this numerical artefact, we compute the WW3 significant wave height anomaly as , where H ctrl s (x) is the significant wave height of the current-free control run.See Appendix B for details.
Figures 3 and 4 show that h s has a wedge-like pattern in the wake of the vortex resulting from wave focussing and defocussing, with h s > 0 mainly for y > 0 and h s < 0 for y < 0. The pattern is not anti-symmetric about y = 0, and positive anomalies are larger than negative anomalies.These characteristics, which indicate a nonlinear response, are increasingly marked as s and U m increase.Specifically, the parameter controls the degree of nonlinearity and hence of asymmetry.We discuss the two limiting regimes γ ≪ 1 and γ ≫ 1 in §5.
There is good overall agreement between WW3 and MA, even though, in the case s = 10, the parameter δ = 0.447 is only marginally small.The pattern is more diffuse for WW3 than for MA, with a less sharply defined wedge and a nonzero h s over a larger proportion of the domain.We attribute the differences to the finiteness of δ (they are more marked for s = 10, δ = 0.447 than for s = 40, δ = 0.224), and to the limited spectral resolution of WW3 (simulations with degraded angular resolution lead to an even more diffuse h s ).The most conspicuous differences between WW3 and MA appear in the scattering region, where the non-zero h s obtained with WW3 appears to contradict the MA prediction that h s = 0.The non-zero h s results from O(ε, δ) terms neglected by MA.Relaxing some of the approximations leading to (3.24) gives a heuristic correction to MA that captures the bulk of the difference with WW3 in the scattering region.We explain this in Appendix C. As further demonstration of the MA approach, we provide a Jupyter notebook accessible at https://www.cambridge.org/S0022112023006869/JFM-Notebooks/files/Figure-3.ipynb,where users can customize the form of the current and the incoming wave spectrum to experiment with the resulting A(x, y, k, θ) dk and h s .

4.2.
Vortex dipole A striking feature of the far-field spectrum and hence of H s is that, according to MA, they depend on the flow only through the deflection ∆(y) in (3.17), proportional to the integral of the vorticity along the direction of dominant wave propagation (the x-direction in our set up).This implies that if the integrated vorticity vanishes because of cancellations between positive and negative contributions, the differences between far-field and incident fields are asymptotically small.This can be interpreted as a form of 'vortex cloaking', whereby an observer positioned well downstream of a flow feature is unable to detect its presence through changes in wave statistics.We demonstrate this phenomenon by examining the scattering of swell by vortex dipoles.
We consider two cases, corresponding to dipoles whose axes (the vector joining the centres of positive and negative vorticity) are, respectively, perpendicular and parallel to the direction of wave propagation.The corresponding vorticity fields are chosen, up to a constant multiple, as the derivative of the Gaussian profile (4.1) with respect to y or x. Figure 5 shows the significant wave height anomaly obtained for the incident spectrum of §4.1 with s = 40 and dipoles with maximum velocity U m = 0.8 m s −1 .
When the dipole axis is in the y-direction (top row) the deflection ∆(y) does not vanish identically.As a result, H s is affected by the flow, strongly for our choice of parameters.This applies to both the MA and WW3 predictions which match closely in the far field.When the dipole axis is in the x-direction (bottom row), ∆(y) = 0.The MA prediction is then that H s = H s⋆ , i.e. h s = 0, everywhere.The WW3 simulation is consistent with this, with only a weak signal in h s .
In general, for a dipole with axis making an angle α with the direction of wave propagation, the deflection ∆(y) is proportional to sin α and the cloaking effect is partial unless α = 0.

Limiting cases
In this section, we return to the far-field asymptotics (3.22) for A in terms of the scaled dependent variables in order to examine two limiting regimes characterized by extreme values of γ = ε/δ.The regime γ ≪ 1 corresponds to a weak flow and/or relatively broad spectrum, leading to a linear dependence of h s on the currents.The opposite regime γ ≫ 1 corresponds to strong flow and/or highly directional spectrum.The wave response is then highly nonlinear in the currents and, as we show below, controlled by the caustics that exist for pure-plane incident waves (γ = ∞).Heller et al. (2008)'s 'freak index', given by ε 2/3 /δ, is the analogue of γ for spatially extended, random currents.
5.1.Linear regime: γ ≪ 1 For γ ≪ 1, we can expand (3.22) in Taylor series to obtain (5.1) This indicates that the flow induces the small correction −γ∆(y−XΘ)∂ Θ A ⋆ (K, Θ) to the action of the incident wave.We deduce an approximation for H s by integrating (5.1) with respect to K and Θ to obtain H 2 s followed by a Taylor expansion of a square root.Alternatively, we can carry out a Taylor expansion of the far-field approximation (3.29) of H s , treating ∆(y) as small.The result is best expressed in terms of the anomaly h s , found to be after reverting to the unscaled variables and setting γ = 1.This simple expression is readily evaluated once the flow, hence ∆(y), and directional spectrum D ⋆ (θ) are specified.For the Gaussian vortex of §4.1 and the directional spectrum in (2.6), the integration can be carried out explicitly, yielding (5.3) This formula makes it plain that h s depends on space through (x/ √ s, y), is antisymmetric about the x axis, and is maximised along the curves y = ± r 2 v + 2x 2 /s.Decay as |x| → ∞ is slowest along these curves and proportional to x −1 .
We illustrate (5.3) and assess its range of validity by comparing it with MA for two sets of parameters in figure 6.The match is very good for s = 10 and U m = 0.4 m s −1 (top row), corresponding to γ = 0.112.It is less good for s = 40 and U m = 0.8 m s −1 , unsurprisingly since γ = 0.447 is not particularly small and the MA prediction is obviously far from linear, with a pronounced asymmetry.The curves y = ± r 2 v + 2x 2 /s shown in the figure are useful indicators of the structure of h s for small enough γ. 5.2.Caustic regime: γ ≫ 1 The limit γ → ∞ corresponds to an incident wave field that is almost a plane wave.It is natural to rescale variables according to Θ → γΘ and X → γ −1 X so that (3.22) becomes A(X, y, K, Θ) = A ⋆ (K, γS(X, y, Θ)) , (5.4) where S(X, y, Θ) (5.5) In (X, y, Θ)-space, the K-integrated action is concentrated in a thin O(γ −1 ) layer around the surface S(X, y, Θ) = 0. Quantities such as H s obtained by further integrating the action with respect to Θ can be obtained by approximating the dependence of right-hand side of (5.4) on S by δ(S).This fails, however, when (X, y, Θ) satisfy both S(X, y, Θ) = 0, and (5.6) The corresponding curves in the (X, y) plane are caustics near which A(X, y, K, Θ) dKdΘ is an order γ 1/2 larger than elsewhere; correspondingly, H s = O(γ 1/4 ).In figure 7 the two caustics meet at a cusp point from opposite sides of a common tangent.The cusp point is located by the condition ∂ 2 Θ S = 0 and the integrated action at the cusp point is O(γ 2/3 ) so that H s = O(γ 1/3 ).We have numerically verified these γ-scalings at the caustics and at the cusp point by varying s in the MA solutions.
For the Gaussian vortex (4.1), the system (5.6) can be solved to obtain an explicit equation for the caustics.This equation is derived in Appendix D and given by (D 6).It describes two curves y(x) emanating from the cusp point at   x = x c given by (D 5).The caustics (which depend on U m but not on s) are indicated on the right panels of figure 3.For the parameters of the figure, the caustics do not map regions of particularly large h s .This is unsurprising since γ is at most 0.447.
To assess how large γ or equivalently s need to be for caustics to be the dominant feature of H s , we show in figure 7 h s computed from MA for U m = 0.8 m s −1 and s = 200 (left panel, γ = 1) and s = 4000 (right panel, γ = 4.47).It is only for s = 4000 that the caustics are evidently controlling the significant wave height pattern.We emphasise that s = 200 and a fortiori s = 4000 are unrealistically large values: observational estimates for s in the open ocean seldom exceed s = 80.We conclude that caustics are unlikely to play a role in real ocean conditions.
With academic rather than practical interest in mind, then, we show in figure 8 the integrated action A dk as a function of y for different three different values of x (identified by dashed vertical lines in figure 7).The figure illustrates how caustics emerge from a fold singularity in the surface S(x, y, θ) = 0 along which action is concentrated in the (x, y, θ) phase space.For x = r v , the surface is a graph over (x, y) and there are no caustics; for x = x c ≈ 3r v , the surface has a single point of vertical tangency (P1 in panel (f) of 7) corresponding to the birth of caustics at a cusp in the (x, y)-plane; for x = 5r c , there are two points of vertical tangency, P2 and P3 in panel (h), corresponding to the two caustic curves.The picture is increasingly blurred as s decreases (compare the right panel of figure 8 with the left panels and with figure 2), explaining the diminishing importance of caustics for H s .

Discussion and conclusion
The main results in this study are obtained by approximate solution of the wave action equation in the four-dimensional position-wavenumber space.A main organizing principle identified by the analysis is that scattering of SGWs by spatially compact currents results in the deflection function, ∆(y) in (3.17).Although ∆ varies linearly with the vertical vorticity of the currents, ∆ figures in a nonlinear transformation of the action density.This nonlinear transformation produces the modulation of the significant wave height H s behind the scattering region, e.g. the expression for H s in (3.29).Quantities that depend on other moments (e.g., Stokes drift) behave similarly and could be readily inferred from our explicit forms (3.15) and (3.22) for the wave action density.
While we have obtained these results for deep-water SGWs, they apply essentially unchanged to other two-dimensional waves with isotropic dispersion relation such as finite-depth SGWs or Poincaré waves.The conclusions we draw about H s can also be rephrased in terms of other root-mean-square quantities relevant to waves other than SGWs.With a little effort, the approach we adopt, based on the matched asymptotics treatment of the wave action equation, could be further extended to three-dimensional waves and to anisotropic dispersion relations.Our results could easily be extended to account for vertically sheared currents using the modified dispersion relation of Kirby & Chen (1989) (which involves a Doppler shift term that is nonlinear in k).
In addition to the WKB approximation used to derive the action conservation equation (2.1) there are two independent approximations involved: (a) the current speed is much less than the group velocity of the incident swell; (b) swell with small directional spreading is incident on a region of spatially compact currents e.g. an axisymmetric vortex or a vortex dipole.Provided that (a) and (b) are satisfied the approximate solution of the wave action equation compares well with numerical solutions provided by WW3.
Approximation (a) is usually justified.To challenge (a) one must consider current speeds such as 2 m s −1 e.g.observed as a peak current speed in the Agulhas system (Quilfen & Chapron 2019).Swell with 100 m wavelength has group velocity ∼ 6 m s −1 so that the small parameter in (a) is as large as 1/3.In less extreme cases approximation (a) will be satisfied.
Approximation (b) is less secure: ocean swell is not sufficiently unidirectional to strongly justify (b) e.g.see the δ-column in table 1.Over long distances, the continuous scattering by uncorrelated currents leads to a broadening of the angular spectrum.When approximation (a) applies, this broadening is described by the directional diffusion equation for wave action derived by Villas Bôas & Young (2020).This diffusion process is one of the mechanisms that makes swell with very small values of δ unlikely.However, our computations for a Gaussian vortex indicate that our asymptotic results are reliable for the moderately small values of δ typical of swell.
Because of the relatively large directional spreading of ocean swell the mathematical ideal of a sharp wave caustic is not realized.Instead the caustic singularity is 'washed out' (Heller et al. 2008).Behind a vortex we find instead an elongated streaky pattern in H s .
Our results show that H s behind an axisymmetric vortex with parameters in table 1 has spatial variation as large as ±30% of the incident constant value H s⋆ .Spatial inhomogeneities in H s of this magnitude are important for wave breaking and exchange of momentum, heat and gas between the ocean and atmosphere.For example, airborne observations of the ocean surface by Romero et al. (2017) indicate that ±30% variations in H s are associated with an order of magnitude increase in whitecap coverage.
The directional diffusion equation of Villas Bôas & Young (2020) uses only approximation (a).One does not need to assume that the wave field is strongly unidirectional or that the currents are spatially compact.Moreover the directional diffusion equation is obtained without detailed consideration of the perturbations to the action spectrum that accompany wave scattering.But there is useful information hiding in these unexamined perturbations to the action spectrum.We are currently engaged in extracting these perturbations, calculating the attendant spatial variability to H s , and relating the statistics of these fluctuations in H s to those of the surface currents.These future developments promise to explain numerical experiments that identify relations between the spectral slopes of surface-current spectra and those of significant wave height (Villas Boâs et al.

2020).
12 s, and minimum source term time step of 5 s.We verified that decreasing the time stepping or the spatial grid spacing does not significantly change the results (not shown).
All simulations are initialized with the narrow-banded wave spectrum in (A 1).Waves enter the domain from the left boundary with initial mean direction θ = 0 • (propagating from left to right), directional spreading parameter s = 10 or s = 10, peak frequency σ ⋆ = 0.61 rad s −1 (peak period of 10.3 s), spectral width δ σ = 0.04, and H s⋆ = 1 m.The boundary condition at the left boundary is kept constant throughout the experiment and each experiment is run until steady state is reached.
As mentioned in §4.1, a control run is conducted in the absence of currents.Although there is no scattering from the currents, a nonuniform h ctrl s = H ctrl s − H s⋆ arises, due to the limited domain size in y, which leads to a reduction of incident wave action from waves arriving from large |y| -an effect that is more pronounced at large x.As s increases, the action density in the incident spectrum is more concentrated in the eastward direction, leading to less leakage of wave action through the top and bottom boundaries and a more spatially uniform h ctrl s .This leakage of wave action corresponds to a reduction of 5% in h ctrl s for s = 10, and 2% for s = 40 towards the right-hand side boundary.

Appendix C. MA-WW3 mismatch in the scattering region
We develop a heuristic correction to MA that we show captures the non-zero h s in the scattering region.First, we note that the non-zero h s in the scattering region from WW3 appears localized, likely caused by the term proportional to ∂ k A in (3.9), as the terms proportional to ∂ θ A result in non-local effects.This observation is confirmed by a WW3 run, which we refer to as WW3 − , where the term in ∂ k A is suppressed in the wave action equation, yielding a more uniform h s in the scattering region (see panel (d) in Figure 9).We then recall that in the MA solution, the insignificance of the ∂ k A term is due to the approximation of a single dominant wavenumber in the steps leading to (3.24).We thus return to the approximation (3.12) of the wave-action transport equation in the scattering region and relax the approximation of replacing k by k ⋆ .We focus on the θintegrated action B(x, k) = A(x, k) dθ.
(C 3) The significant wave height is deduced by integration as  We emphasise the heuristic nature of this approximation (MA + ) which is formally no more accurate than the MA approximation H s (x) = H s⋆ since it neglects some, though not all, O(δ) terms.Nonetheless, it captures most of the significant wave height anomaly close to the Gaussian vortex, as figure 9 demonstrates under parameters s = 40 and U m = 0.8 m s −1 .

Figure 2 :
Figure 2: Wavenumber-integrated action density A(x, y, k, θ) dk as a function of y and θ at x = −5 rv , 0, rv , 3 rv and 5 rv from WW3 (left) and MA (Eq.(4.5), right) for swell impinging on a Gaussian vortex with Um = 0.8 m s −1 .The directional spreading of the incident spectrum is s = 40.

Figure 3 :Figure 4 :
Figure 3: Significant wave height anomaly hs(x, y) from WW3 (left column) and MA (right column) for swell impinging on a Gaussian vortex.Each row corresponds to the indicated values of the directional spreading parameter s of the incident wave spectrum and of the maximum velocity Um (in m s −1 ).The corresponding non-dimensional parameters are given in Table 1.The dashed circles has radius rv around vortex center.The solid lines on the right panels indicate the caustics computed from (D 6).The colourbars differ between rows but are the same within each row.White corresponds to hs = 0 in all panels.The customizable notebook that generates panel (h) by default can be accessed at https://www.cambridge.org/S0022112023006869/JFM-Notebooks/files/Figure-3.ipynb.

Figure 5 :
Figure 5: Swell impinging on vortex dipoles with axes perpendicular (top) and parallel (bottom) to the dominant direction of wave propagation (x-axis).The vorticity (colour) and velocity (vectors) are shown (left) together with the significant wave height anomaly hs from WW3 (middle) and MA (right).The directional spreading parameter s = 40 and the maximum flow velocity is 0.8 m s −1 .

Figure 6 :
Figure 6: Significant wave height anomaly hs(x, y) for swell impinging on a Gaussian vortex: comparison between the predictions of MA (left) and and its γ → 0 limit ((5.3), right column).The set up is as in figure 3 with parameters s and Um (in m s −1 ) as indicated.Dashed lines indicates the curves y = ± r 2 v + 2x 2 /s where hs reach maximum amplitudes according to (5.3).

Figure 7 :
Figure 7: Caustics for swell impinging on a Gaussian vortex: the caustics (D 6) (solid lines) are superimposed to the MA prediction of hs for Um = 0.8 m s −1 and the indicated values of s.The dashed vertical lines correspond to the values of x = rv, 3 rv and 5 rv used in figure 8.

Figure 8 :
Figure 8: Wavenumber-integrated action density A(x, y, k, θ)dk as a function of y and θ for x = rv, 3 rv and 5 rv corresponding to the significant wave height shown in figure 7 for s = 200 (left column) and s = 4000 (right column).P1 in panel d corresponds to the values of (x, y) of the cusp from where the caustics emanate; P2 and P3 are associated to points on each of the two caustics.
4)We now change the integration variable, taking advantage of the localisation of B ⋆ (k) to ignore the corresponding change in the lower limit of integration and

Figure 9 :
Figure 9: Significant wave height anomaly hs computed from (a) WW3, (b) MA (same as figures 3g,h), (d) WW3 − (where the term proportional to ∂ k A is switched off) and (e) MA + as in (C6).Panels (c) and (f) show the differences between (a) and (b), and (a) and (e), respectively.All plots have the same colour bar.

Table 1 :
.1) obtained with s Um (m s −1 ) δ = 2/s ε = Um/c⋆ γ = ε/δ Parameters corresponding to each configuration in section 4.1, arranged in the order of the rows in figure3.In all cases the group speed is c⋆ = 8 m s −1 , corresponding to a 166 m wavelength and 10.3 s period.Um is the maximum vortex velocity and the vortex radius is rv = 25 km.theWave Height, Water Depth, and Current Hindcasting third generation wave model (WAVEWATCH III, hereafter WW3).The incident spectrum used for WW3 is described in Appendix A. The directional function for this spectrum is the Longuet-Higgins et al. (1963) model