Interaction of linear modulated waves with unsteady dispersive hydrodynamic states

A new type of wave-mean flow interaction is identified and studied in which a small-amplitude, linear, dispersive modulated wave propagates through an evolving, nonlinear, large-scale fluid state such as an expansion (rarefaction) wave or a dispersive shock wave (undular bore). The Korteweg-de Vries (KdV) equation is considered as a prototypical example of dynamic wavepacket-mean flow interaction. Modulation equations are derived for the coupling between linear wave modulations and a nonlinear mean flow. These equations admit a particular class of solutions that describe the transmission or trapping of a linear wave packet by an unsteady hydrodynamic state. Two adiabatic invariants of motion are identified that determine the transmission, trapping conditions and show that wavepackets incident upon smooth expansion waves or compressive, rapidly oscillating dispersive shock waves exhibit so-called hydrodynamic reciprocity recently described in Phys. Rev. Lett. 120, 144101 (2018) in the context of hydrodynamic soliton tunnelling. The modulation theory results are in excellent agreement with direct numerical simulations of full KdV dynamics. The integrability of the KdV equation is not invoked so these results can be extended to other nonlinear dispersive fluid mechanic models.


Introduction
The interaction of waves with a mean flow is a fundamental and longstanding problem of fluid mechanics with numerous applications in geophysical fluids (see e.g. Mei et al. (2005), Bühler (2009) and references therein). The key to the study of such an interaction is the scale separation, whereby the length and time scales of the waves are much shorter than those of the mean flow. For linear waves, the mean flow is not affected by the amplitude variations (i.e. the effects of induced mean flow are neglected) and can be specified in advance (see, e.g., Peregrine (1976)). Once the dynamical equations are linearised, the variations of the mean flow enter the coefficients of the linear system, which is formally equivalent to the dynamics of linear waves in non-uniform and unsteady media.
Due to the multi-scale character of wave-mean flow interaction, a natural mathematical framework for its description is Whitham modulation theory Whitham (1965aWhitham ( , 1999. Although the initial motivation behind modulation theory was the study of finite-amplitude waves, it was recognised that the wave action equation that plays a fundamental role in Whitham theory (Hayes, 1970) was also useful for the study of linearised waves on a mean flow, see e.g. Garrett (1968) and Grimshaw (1984). It was used in Bretherton & Garrett (1968) and Bretherton (1968) to examine the interaction between short-scale, small amplitude internal waves and a mean flow in inhomogeneous, moving media. The outcome of this pioneering work was the determination of the variations of the wavenumber, frequency and amplitude of the linearised wavetrain along group velocity lines. Subsequently, this work was extended in Grimshaw (1975) to finite amplitude waves, incorporating the perturbative effects of friction and compressibility, as well as the leading order effect of rotation.
The modulation theory of linear wavetrains in weakly non-homogeneous and weakly non-stationary media (where "weakly" is understood as slowly varying in time and/or space) was developed in Whitham (1965a) and Bretherton & Garrett (1968). It was shown that the modulation system for the wavenumber k, the frequency ω and the amplitude a is generically composed of the conservation equations where the dispersion relation ω(k; α(x, t)) and the wave action density A(a, k; α(x, t)) depend on the system under study with α(x, t) being a set of slowly varying coefficients describing nonhomogeneous non-stationary media, including the effects of the prescribed mean flow, e.g., the current. Equations (1) were applied to the description of the interaction of water waves with a steady current in (Longuet-Higgins & Stewart, 1961;Peregrine, 1976;Phillips, 1980;Peregrine & Jonsson, 1983;Whitham, 1999;Mei et al., 2005;Bühler, 2009). We briefly outline some classical results from the above references relevant to the developments in this paper. Consider a deep water surface wave interacting with a given non-uniform but steady current profile U (x). Assuming slow dependence of U on x, the linear dispersion relation reads ω = U (x)k + σ(k) where σ(k) = √ gk is the so-called intrinsic frequency, i.e. the frequency of the wave in the reference frame moving with the current, and the wave action density has the form A = ρga 2 /σ(k). Since U only depends on x, we look for a steady solution of the modulation equations (1), which yield ω x = 0 and (∂ k ωA) x = 0.
Suppose further that U (x) slowly varies between U − = 0 and U + < 0. The wavenumber of the linear wave then slowly changes from some k − to k + , and the conservation of the frequency yields the following relation: gk − = gk + + U + k + , thus the linear wave shortens, k + > k − , when it propagates against the current (U + < 0). Interestingly, if k − = g/16U 2 + , the group velocity ∂ k ω vanishes for k = k + = g/4U 2 + , and no energy can propagate against the current, i.e. the wave is "stopped", or "blocked" by the current (Taylor, 1955;Lai et al., 1989). The amplitude of the linear wave is given by the conservation of the wave action flux: a 2 ∝ 1/(U +∂ k σ). Since U + < 0, the wave amplitude increases a + > a − . If k − is close to g/16U 2 + , the amplitude becomes extremely large and the wave breaks. As a matter of fact, the linear approximation fails to be valid for such waves. As noted in Peregrine (1976), "such a stopping velocity ... leads to very rough water surfaces as the wave energy density increases substantially. Upstream of such points, especially if the current slackens, the surface of the water is especially smooth as all short waves are eliminated." This phenomenon has been observed when the sea draws back at the ebb of the tide where an opposing current increases wave steepness and, as a result, wave breaking occurs (see for instance Johnson (1947)). It also enters the mechanism of some pneumatic and hydraulic breakwater (cf., e.g., Evans (1955)) where the injection of a local current destabilises the waves and prevents them from reaching the shore. More recently wave blocking has been used to engineer the so-called white hole horizon for surface waves in the context of analogue gravity (Rousseaux et al., 2008(Rousseaux et al., , 2010. Similar problems for a non-uniform and unsteady mean flow have also been studied. The necessity to consider the media unsteadiness has been first recognised in Unna (1941); Barber (1949); Longuet-Higgins & Stewart (1960). Due to the non-stationary character of the problem, the frequency, as well as the wavenumber, are not constant. Various "unsteady set ups" have been studied using the linear theory (1), such as the influence of an unsteady gravity constant on water waves (Irvine, 1985), the effect of internal waves on surface waves (Hughes, 1978), the water wavetidal wave interaction (Tolman, 1990), or the influence of current standing waves on water waves (Haller & TubaÖzkan-Haller, 2007).
In all described examples, the mean flow or the medium nonhomogeneity were prescribed externally, which results in a simple modulation system (1), consisting of just two equations with variable coefficients, with several straightforward implications for the wave's wavelength and amplitude as outlined above. In this work, we study a different kind of wave-mean flow interaction, where the mean flow dynamically evolves in space-time so that the variations of both the wavetrain and the mean flow are governed by the same nonlinear dispersive PDE but occur in differing amplitudefrequency domains. The dynamics of the small amplitude short-wavelength wave are dominated by dispersive effects while the large-scale mean flow variations are dominated by nonlinearity. In this scenario, the modulation system (1) for the linear wave couples to an extra nonlinear evolution equation for the mean flow. The form of the mean flow equation depends on the nature of the largescale unsteady fluid state involved in the interaction. For the simplest case of a smooth expansion (rarefaction) wave, the mean flow equation coincides with the long-wave, dispersionless, limit of the original dispersive PDE. However, if the large-scale, nonlinear state is oscillatory, as happens in a dispersive shock wave (undular bore), the derivation of the mean flow equation requires full nonlinear modulation analysis originally presented in Gurevich & Pitaevskii (1974) (see also  and references therein). We show that in both cases, the wave-mean flow interaction exhibits two adiabatic invariants of motion that govern the variations of the wavenumber and the amplitude in the linear wavetrain, and prescribe its transmission or trapping inside the hydrodynamic state: either a rarefaction wave (RW) or a dispersive shock wave (DSW). Trapping generalises the aforementioned discussion of blocking phenomena to time-dependent, nonlinear mean flows.
As a basic prototypical example, we consider dynamic wavepacket-mean flow interactions in the framework of the KdV equation The settings we consider are illustrated in Fig. 1. The linear wavepacket propagating with group velocity −3k 2 relative to the background, say u = u 0 , is incident from the right upon an unsteady dispersive-hydrodynamic state: a RW or a DSW. We derive modulation equations describing the coupling between the modulations in the linear wave packet and the variations of the background mean flow and show that the linear wave is either transmitted through or trapped inside the unsteady hydrodynamic state. The transmission/trapping conditions are determined by two adiabatic invariants of motion that coincide with Riemann invariants of the modulation system on a certain integral surface. The paper is organised as follows. In Sec. 2 we introduce the mean field approximation and use the Whitham averaged Lagrangian approach to derive a general modulation system describing the interaction of a linear modulated wave with a nonlinear dispersive hydrodynamic state: either RW or DSW. This system consists of two usual modulation equations describing conservation of wave number and of the wave action, which are coupled to the simple wave evolution equation describing the mean flow variations in the RW/DSW. The characteristic velocity of the mean flow equation depends on the nature -smooth (RW) or oscillatory (DSW) -of the hydrodynamic state. The obtained full modulation system, despite being nonstrictly hyperbolic, is shown to generically possess a Riemann invariant associated with the linear group velocity characteristic, in addition to the simple wave Riemann invariant associated with the mean flow evolution. Moreover, we show that the wave action modulation equation can be further manipulated in the diagonal form effectively exhibiting an extra Riemann invariant on a certain integral surface.
In Sec. 3, we consider the model problem of the plane wave-mean flow interaction whereby the mean flow variations are initiated by Riemann step initial data. Within this framework, the Riemann invariants of the modulation system found in Sec. 2 are shown to play the role of adiabatic invariants of motion which predict the transmission conditions through the RW. The transmission through DSW is then determined by the same conditions as in the RW case following the socalled hydrodynamic reciprocity recently described in the context of soliton-mean flow interactions (Maiden et al., 2018).
The results of Sec. 3 are employed in Sec. 4 to study the physically relevant-and numerically tractable-case of the interaction of localised wavepackets with RWs and DSWs. To study the wavepacket-mean flow interaction, we introduce a partial Riemann problem (wavepacket-RW interaction) and show that the variation of the wavepacket dominant wavenumber is governed by the conservation of the adiabatic invariant identified in Sec. 3 and thus yields the same transmission and trapping conditions for the wave packet as in the full Riemann problem (plane wave-RW interaction). The same conditions are valid, via hydrodynamic reciprocity, for the wavepacket-DSW interaction case. Wavepacket trajectories inside the RW and the DSW are also determined analytically and compared with the results of numerical resolution of the corresponding partial Riemann problem. It is shown that wavepackets exhibit a speed shift, as well as a phase shift, as a result of their interaction with a hydrodynamic state.
In Sec. 5, we draw conclusions and identify applications and perspectives of further development of our work. Appendix A describes the numerical implementation of the partial Riemann problem employed in Sec. 4.
2 Modulation dynamics of the linear wave-mean flow interaction 2.1 Mean field approximation and the modulation equations In this section, we shall introduce the mean field approximation that enables a straightforward derivation of the modulation system describing linear wave-mean flow interaction. The full justification of this approximation for the case of the interaction with a RW can be done in the framework of standard multiple-scales analysis (Luke, 1966), equivalent to single-phase modulation theory (Whitham, 1999). The justification for linear wave-DSW interaction is more subtle, requiring the derivation of multiphase (two-phase) nonlinear modulation equations Ablowitz & Benney (1970) and making linearisation in one of the oscillatory phases. To avoid unnecessary technicalities, we simply postulate the approximation used and then justify its validity by comparison of the obtained results with direct numerical simulations of the KdV equation.
To describe the interaction of a linear dispersive wave with an extended nonlinear dispersivehydrodynamic state (RW or DSW), we represent the solution u(x, t) of the KdV equation (2) as a superposition where u H.S. (x, t) corresponds to the RW or DSW solution, and ϕ(x, t) corresponds to a small amplitude field describing the linear wave. In order to extract the dynamics of ϕ(x, t), we make the mean field (scale separation) approximation by assuming that u H.S. (x, t) is locally (i.e. on the scale ∆x ∼ ∆t = O(1)) periodic and replace the dispersive hydrodynamic wave field u H.S. (x, t) with its local mean (period average) value u(x, t) while assuming that ϕ u, u x /u ϕ x /ϕ . For a smooth, slowly varying hydrodynamic state (RW) such a replacement is natural since locally one has u H.S. = u, but for the oscillatory solutions describing slowly modulated nonlinear wavetrains in a DSW, u H.S. (x, t) = u so the mean field approximation would require justification via a careful multiple scale analysis (Ablowitz & Benney, 1970). In particular, a detailed analysis of possible resonances between the DSW and the wavepacket will be necessary (Dobrokhotov & Maslov, 1981). Such a mathematical justification will be the subject of a separate work, while here we shall postulate the outlined mean field approximation and show that it enables a remarkably accurate description of the linear field ϕ, which can be thought of as propagating on top of the mean flow.
Within the proposed mean field approximation, the small amplitude wave field ϕ satisfies the linearised, variable coefficient KdV equation where for the case of linear wave-RW interaction, the mean flow u(x, t) evolves according to the Hopf equation For the propagation of a linear wave through a DSW, which is a rapidly oscillating nonlinear structure, the linearised field ϕ satisfies the same equation (4) but the mean flow evolution is no longer described by the Hopf equation (5). Instead, its dynamics are approximately governed by the modulation equation (Gurevich & Pitaevskii, 1974)  with V (u) given parametrically by where K(m) and E(m) are the complete elliptic integrals of the first and the second kind respectively (Abramowitz & Stegun, 1972), m ∈ [0, 1], and u − > u + are the values of u at the left and right constant states respectively, connected by the DSW; Fig. 2 displays the variation of the characteristic speed V (u) for (u − , u + ) = (1, 0). The mean flow variation u(x, t) in the DSW is determined by the self-similar solution V (u) = x/t and is displayed in Fig. 2. Equations (6), (7) follow from the Whitham modulation system obtained by averaging of the KdV equation over the family of nonlinear periodic (cnoidal wave) KdV solutions (Whitham, 1965b). This system consists of three hyperbolic equations that can be diagonalised in Riemann invariant form. The DSW modulation represents a simple wave (more specifically, a 2-wave, see e.g. ) solution of the Whitham equations, in which two of the Riemann invariants are set constant to provide continuous matching with the external constant states u ± (Gurevich & Pitaevskii, 1974), see also (Kamchatnov, 2000). As we shall show, Eqs. (4), (6), (7) provide an accurate description of the interaction between the linear wave and the DSW so that the dynamics of ϕ(x, t) are predominantly governed by the variations of the DSW mean value u(x, t). In what follows, we will be using the general equation (6) for the mean flow, assuming that V (u) = u for linear wave-RW interaction. We note that a similar mean flow approach, in which the oscillatory DSW field was replaced by its mean u has been recently successfully applied to the description of a soliton-DSW interaction in Maiden et al. (2018).
Equations (4), (6) form our basic mathematical model for linear wave-mean flow interaction. We shall proceed by constructing modulation equations for this system. One may question the wisdom of incorporating the decoupled mean flow equation (6) to Eq. (4) rather than simply prescribing an arbitrary mean flow externally. As we will see, the mathematical structure of Eqs. (4), (6) enables a convenient solution that is not available for generic mean flows u(x, t). Moreover, Eqs (4), (6) transparently reveal the multiscale structure of the dynamics: a "fast" equation (4) for the linear waves and a "slow" equation (6) for the mean flow.
Let ϕ(x, t) describe a slowly varying wavepacket: where we assume: a x /a, k x /k, ω x /ω k , a t /a, k t /k, ω t /ω ω .
Assuming the slowly varying background u, i.e. u x /u ϕ x /ϕ, we determine the modulation equations for a and k through Whitham's variational approach (Whitham, 1999). The Lagrangian associated with Eq. (4) is The average of the Lagrangian L over one oscillation period 2π leads to the averaged Lagrangian The variation of the mean field u is assumed to be negligible over a single wave period, θ ∈ [0; 2π], but it includes an explicit slow dependence on x, t in the averaged Lagrangian (11). The modulation equations for the linear wave interacting with the hydrodynamic state are then equivalent to the modulation equations for a wavepacket propagating through a non-uniform medium (Whitham, 1999) and are obtained from the variational principle for the fields a(x, t) and θ(x, t). It is convenient to represent the averaged Lagrangian (11) in a general form with α(x, t) being a vector of parameters depending explicitly on x and t; in the present case (11), we have α = u(x, t) a scalar parameter and G = − ω k + α − k 2 . The consistency condition θ xt = θ tx yields the wave conservation equation while the remaining two equations are obtained by considering variations in δa and δθ to give the standard modulation system δa : δθ : which remain valid in the presence of explicit x, t-dependence in the average Lagrangian L .
which is recognised as the linear dispersion relation ω = Ω(k; α) in the nonhomogeneous medium; for the KdV equation (2), Ω(k; u) = u k − k 3 . Taking into account (13), the wave action equation (16) becomes where v g = Ω k is the group velocity and A = G ω a 2 is the wave action (cf. Eq. (1)). For the KdV equation, we obtain from (18), (11) ∂ ∂t which can be cast in a physically transparent amplitude form by expanding (19) and using (14) ( with the last term describing the wave-mean flow interaction. The modulation equations (6), (14) and (20) are the main subject of this paper. For V (u) = u in eq. (6), they can be obtained from the single-phase KdV-Whitham system for modulated cnoidal waves in the small-amplitude, zero induced mean flow approximation, i.e., by neglecting the term proportional to a 2 (a 2 ) x in the amplitude equation and the term proportional to (a 2 ) x in the equation for u and (see Whitham (1999) Ch.16). However, as we already mentioned, the derivation of this system for the interaction of a wavepacket with a DSW, where V (u) is given by Eq. (7), is a non-trivial, technical task requiring complicated multiple scale analysis Ablowitz & Benney (1970) or multi-phase Whitham theory Flaschka et al. (1980). The formal direct mean field Lagrangian dynamics approach adopted here is more convenient for our purposes and also admits broad generalisations. We shall justify its validity by making comparisons with direct numerical simulations. Finally, we derive a useful consequence of the wave conservation law (14) for a "wavepacket train" consisting of a superposition of two slowly modulated plane waves with close wavenumbers k and k + δk where δk k, which corresponds to beating of the two waves. The conservation of waves for these two waves read Hence, for δk k, the subtraction of these two equations reduces to a conservation equation very similar to (19): Concluding this section, we note that the modulation equations (14), (19) are quite simple and definitively not new. However, unlike in previous studies, they are now coupled to the mean field equation (6). As we shall show, the system consisting of Eqs. (6), (14), (19), and (22) equipped with appropriate initial conditions, yields straightforward yet highly non-trivial implications, especially in the case of wavepacket-DSW interaction, which is very difficult to tackle using direct (nonmodulation) analysis.

Riemann invariants
In what follows, we shall use a more natural notation ω(k, u) for the dispersion relation instead of Ω(k; u). Additionally, we shall use an "abstract" field A(x, t) instead of a(x, t) 2 /k(x, t) and δk(x, t) such that the reduced modulation system composed of Eqs. (6), (14) and (19) (or (22)) can be cast in the general form We note that the system (23) has the double characteristic velocity v g and thus is not strictly hyperbolic. The first two equations, (23a) and (23b), are decoupled and can always be diagonalised such that Eq. (23b) assumes the form for the Riemann invariant q = Q(k, u). Generally, the Riemann invariant q as a function of u and k is found by integrating the characteristic differential form For the case of linear wave-RW interaction, we have V (u) = u, so Ξ = −3k 2 dk + kdu which can be integrated after multiplying by the integrating factor 1/k to yield explicit expressions for the Riemann invariant q and the associated characteristic velocity It follows from (24) that q = const ≡ q 0 along the double characteristic dx/dt = v g , which enables one to manipulate equation (23c) into the form valid along dx dt = v g , where with u 0 a constant of integration. The quantity pA thus can be viewed as a Riemann invariant of the system (23) on the integral surface Q(k, u) = q 0 which will prove useful in the analysis that follows. We stress, however, that pA is not a Riemann invariant in the conventional sense since the system (23) is not strictly hyperbolic. See Maiden et al. (2018) for a similar construction in the context of soliton-mean flow interaction. For V (u) = u, the integral in (28) is readily evaluated to give, taking into account (26), Note, that the described diagonalisation of the reduced modulation system (23) for the linear wave-mean flow interaction, unlike the existence of Riemann invariants of the general Whitham system for modulated cnoidal waves (Whitham, 1965b), does not rely on integrability of the KdV equation. In fact, the possibility of this diagonalisation is general and is a direct consequence of the Figure 3: Schematic of the initial conditions for the interaction between a plane wave and a hydrodynamic state.
absence of induced mean flow for linearised waves, so that the dynamics of the wave parameters (k, ω, a) are decoupled from the dynamics of the mean flow u. Here, we reap the benefits of jointly considering the evolution of the mean flow, wavenumber conservation, and the field equation in (23) by recognising that they can be cast in diagonal, Riemann invariant form (23a), (24), and (27) along 3 Plane wave-mean flow interaction: the generalised Riemann problem

Adiabatic invariants and transmission conditions
Before studying the interaction of localised wavepackets with mean flow in the framework of the basic system (4), (6) we consider a model problem of the unidirectional "scattering" of a linear plane wave (PW) by a nonlinear a hydrodynamic state (RW or DSW) initiated by a step in u. We denote the incident PW parameters at x → +∞ as k + , a + and the transmitted PW parameters at x → −∞ as k − , a − . To find the transmission relations, we consider the generalised Riemann problem (see Fig. 3) We call this Riemann problem "generalised" as it is formulated for the modulation system rather than for the original dispersive model (4), (6).
In the interaction of a linear PW with both a RW and a DSW, the evolution of u(x, t) is described by the self-similar expansion fan solution of the mean flow equation (6): while the Riemann invariants q and pA are constant throughout, implying the relations where A ± = a 2 ± /k ± . The expressions for Q(k, u) and P (k, u) in the PW-RW interactions for KdV are given by Eqs. (26) and (29), respectively. For PW-DSW interaction, when V (u) is given by Eq. (7), simple explicit expressions for Q and P in terms of k and u are not available. However, they can be obtained by integrating (25) and evaluating (28), e.g., numerically. The edge speeds of the expansion fan (31) are given by RW : DSW : Expressions (35) follow from Eq. (7) upon taking m → 0 + and m → 1 − for the trailing and the leading DSW edge respectively, see Gurevich & Pitaevskii (1974). Given u(x, t) described by (31), the conservation of q and pA in (32), (33) yields not only the PW transmission relations but also the slow variations of the PW parameters k(x/t) and a(x/t) due to the interaction with the mean flow in the hydrodynamic state. Constant q and pA can thus be seen as adiabatic invariants of the PW-mean flow interaction.
The conservation relations (32), (33) also describe wave-mean flow interaction for a system of two beating, superposed slowly modulated plane waves with close wavenumbers k and k + δk interacting with a RW (cf. Sec. 2.1). In this case, the adiabatic variation of k and δk are described respectively by (32) and (33) where A corresponds now to δk.

Plane wave-rarefaction wave interaction
A RW is generated when u − < u + , and the resulting mean flow variation is described by Eq. (31) with characteristic velocity V (u) = u. In this case, explicit expressions for the adiabatic invariants q and pA can be obtained using Eqs. (26) and (29). The conservation relations (32), (33) then yield where the second condition was obtained by using A ± = a 2 ± /k ± . It is remarkable that the interaction of a PW with a non-uniform, unsteady hydrodynamic state does not change the PW amplitude, which is in sharp contrast with the case of water wave-steady mean flow interaction outlined in Sec. 1, where the amplitude varies following the inhomogeneities of the current. Fig. 4 displays the comparison between the relation (36) and the wavenumbers obtained in the numerical simulations of the linear wave-mean flow interaction. In numerical simulations we employed the partial Riemann problem (see Sec. 4.1), where the PW, which is a mathematical abstraction, was replaced by a broad but localised linear wavepacket so that only two variables (u and k) experience an initial jump. This replacement is legitimate here since, as we show in Sec. 4.1, relation (36) is also valid in the localised wavepacket-RW interaction with k + and k − being the dominant wavenumbers of the incident and transmitted wavepackets respectively.
One can see that Eq. (36) yields the transmission condition: the transmitted PW exists if its wavenumber k − is a real number. This requirement reduces to the following condition for the wavenumber k + of the incident PW: The PW wavenumber variations during its interaction with the RW follow from Eqs. (31), (26) and (32): It follows from (39) that the PW wavelength increases inside the RW (since u + − x/t > 0), so if condition (38) is not fulfilled, one has k = 0 for the position i.e., the PW gets absorbed by the RW and does not get transmitted through. This result will receive further interpretation in Sec. 4 in the context of the interaction of a localised wavepacket with a RW as "wave trapping" inside the hydrodynamic state. In that case, a wavepacket would asymptotically, as t → ∞, propagate with the same group velocity as the local RW background, v g (k; u) = x/t = u, which implies that its dominant wavenumber will decay with time.
One can also consider the interaction between two beating superposed PWs and a RW (see the discussion in Sec. 2.1) where the conservation of the adiabatic invariant p(k, u)A with A = δk yields: As was mentioned in the previous section, the beating pattern created by the superposition of two PWs with close wavenumbers (δk k) can be seen as a "wavepacket train" of period L = 4π/δk. Thus Eq. (40) provides the relation between the period of the incident train L + and the transmitted train L − . The difference can be interpreted as the phase shift between the incident and the transmitted wavepackets. Similarly, one can interpret L ± = 4π/δk ± as the widths of the wavepackets before and after transmission. Their relation is then given by k − /L− = k + /L + .

Plane wave-DSW interaction: hydrodynamic reciprocity
We now consider the initial condition (30) with u − > u + resolving into a DSW. In this case, the modulation of the mean flow is described by the simple wave equation (31), (7) and the expressions for the adiabatic invariants q and pA differ from Eqs. (26), (29) obtained for PW-RW interaction. As a result, the conditions (32), (33) for the conservation of q and pA describe a very different adiabatic evolution of the PW parameters inside the dispersive hydrodynamic state. We shall consider this evolution later in Sec. 4.4, while here we describe a very general property of PW interaction with dispersive hydrodynamic states termed hydrodynamic reciprocity, that was initially formulated in Maiden et al. (2018) in the context of the mean field interaction of solitons with dispersive hydrodynamic states.  (1, 0). The markers correspond to the relation determined numerically for the appropriate partial Riemann problem (cf. Sec. 4.1); more specifically, the crosses ( ) are determined when the linear wave interacts with a RW and the circles ( ) when the linear wave interacts with the DSW. The relation between k − and k + is independent of the nature of the hydrodynamic state, as illustrated by the continuous modulation solution used in the construction of the contour plot in the left panel.
We notice that the PW-RW and PW-DSW interactions in the mean flow approximation are described by the solutions of the same Riemann problem considered in t > 0 and t < 0 half-planes. Indeed, the Riemann problem (23), (30) with u − > u + describes PW-DSW interaction for t > 0. However, when considered for t < 0, this problem describes PW-RW interaction, characterised by the transition relations (36) and (37) which also hold for t = 0. Then, continuity of the simple wave modulation solution for all (x, t) (illustrated by Fig. 4.a), except at the origin (x, t) = (0, 0), implies that the transition relations (36) and (37) derived for the PW-RW interaction (t < 0) must also hold for t > 0, i.e., for the PW-DSW interaction. This hydrodynamic reciprocity can be verified in Fig. 4.b where we compare the relations between k − and k + obtained numerically for the evolution of PW-RW and PW-DSW interactions in the full KdV equation. We see the excellent agreement, confirming that relation (36), as well as the transmission condition (38), indeed hold for the PW interaction with both nonlinear dispersive hydrodynamic states: RW and DSW.
The perfect agreement between relation (36) and the numerical solution of the Riemann problem in Fig. 4 also confirms the mean field hypothesis which is behind the basic mathematical model (4) of this paper. Although the mean field assumption is clear for the PW-RW interaction where the hydrodynamic state solution u H.S. , up to small dispersive corrections at the RW corners, coincides with the solution of the Hopf equation (5) for mean flow evolution u, it is no longer the case for the highly nontrivial PW-DSW interaction where u H.S. describes a rapidly oscillating structure, which is radically different from its slowly varying mean flow u satisfying the equations (6), (7).
4 Interaction of linear wavepackets with unsteady hydrodynamic states

Partial Riemann problem
Having considered the model case of PW interaction with dispersive hydrodynamic states, we now proceed with a more physically relevant example of a similar interaction involving localised linear wavepackets instead of PWs. To model such an interaction, the Riemann problem (23), (30) should be modified to take into account the localised nature of the wavepacket. To this end, we introduce the partial Riemann problem where f (x − X 0 ) is a localised amplitude distribution, centered at X 0 and decaying as |x| → ∞. In what follows, the position of the wavepacket, defined as a group velocity line, is denoted by X(t), so that according to (42), X(0) = X 0 , see Eq. (47) below. The quantities k + and k − here denote the dominant wavenumbers of the wavepackets for x > 0 and x < 0 respectively. Note that, although the wavepacket dominant wavenumber is defined only along the group velocity line, we treat it here as a spatio-temporal field k(x, t). See Maiden et al. (2018) for a similar extension made to define a soliton amplitude field in the context of the soliton-mean flow interaction problem. While the formulation (42) assumes the simultaneous presence of two wavepackets at t = 0, only one of them -we shall call it the incident wavepacket -is physically realised due to the localised nature of the amplitude distribution f (X). The additional, fictitious wavepackets will yield all the transmission, trapping information of the incident wavepacket.
The wavepacket is assumed to be sufficiently broad such that its amplitude a(x, t) does not vary significantly over one period 2π/k of the oscillation of the carrier wave, see conditions (9). At the same time, the typical width of the wavepacket is assumed to be much less than the width of the (expanding) hydrodynamic state in order for the mean field approximation to be applicable. This latter condition is readily satisfied for sufficiently large t as the dispersive broadening rate of a linear quasi-monochromatic wavepacket is determined by the difference |v g (k max ) − v g (k min )| ∼ δk, where δk = k max − k min 1 is the spectral width (Whitham, 1999), and is much smaller than the broadening rate of the RW or DSW, which is proportional to |u − − u + | ∼ 1.
The partial Riemann problem (23), (42) implies two possible interaction scenarios: (i) a "rightincident" interaction where a wavepacket, initially placed at X 0 > 0, propagates with group velocity v + g = u + − 3k 2 + and enters either expanding hydrodynamic structure whose leading edge velocity is V (u + ) > v + g (see (34), (35)); (ii) a "left-incident" interaction, where the wavepacket is initially placed at X 0 < 0 so that the interaction only occurs if the velocity of the trailing edge of the hydrodynamic structure is less than the group velocity of the wavepacket, i.e. V (u − ) < v − g = u − − 3k 2 − . It follows from (34), (35) that this can happen only for a DSW but not for a RW. The subsystem (23a), (23b) and (42) for u and k has already been solved in the previous section. The simple wave solution of this problem is given by Eqs. (31) and (32), and thus the relation between k − and k + (36) obtained for PWs, approximately holds for the wavepacket-mean flow interaction. As a result, the wavepacket is subject to the transmission condition (38), and Hydro. state Wavepacket X 0 < 0 • trapped in the DSW otherwise • always transmitted Table 1: Classification of possible configurations for the wavepacket -hydrodynamic state interaction with k ± > 0.
the possible interaction configurations are summarised in Table 1. We now perform numerical verification of the relation (36) by numerically solving the partial Riemann problem (23), (42). In fact, it was the partial Riemann problem that was numerically solved to verify the hydrodynamic reciprocity property illustrated in Fig. 4. The numerical implementation of the partial Riemann problem (23), (42) is detailed in Appendix A.

Conservation of the integral of wave action
We now proceed with the determination of the wavepacket amplitude variation resulting from the interaction with dispersive hydrodynamic states. It is well known that KdV dispersion leads to wavepacket broadening so that its amplitude a(x, t) is not constant but decreases during propagation on a constant mean flow u − or u + in order to conserve the integral ∞ −∞ a 2 dx, which leads to the standard dispersive decay estimate a ∼ t −1/2 for t 1 (Whitham, 1999) . Thus, we cannot expect the amplitude transmission relation (37) derived for PWs to remain valid for localised wavepackets. To address this issue, instead of considering the amplitude of the wave, we consider the integral of wave action between two group lines dX 1,2 /dt = [v g (k(x, t), u(x, t))]| x=X1,2(t) . (43) is conserved during linear wavepacket propagation through a hydrodynamic state with varying mean flow u(x, t) (cf., the conservation of wave action in Whitham (1965a); Bretherton & Garrett (1968)). If the incident wavepacket is transmitted, we evaluate the integral (43) before and after the interaction, i.e. at times t = t + and t = t − when the wave is localized at the right and at the left of the hydrodynamic state, respectively, and where the dominant wavenumber field is uniform k(x, t ± ) = k ± for a(x, t ± ) = 0. Thus, the conservation of the integral of wave action E(t − ) = E(t + ) yields: where we replace the limits of integration X 1 and X 2 by −∞ and +∞ since the wavepacket is localised in space. The relation (44) is valid in both linear wavepacket-RW and -DSW interaction, as illustrated in Fig. 5.

It then follows from (19) that the integral
In the case of a broad wavepacket of almost constant amplitude, we have the following approximation: where a ± and L ± are, respectively, the constant amplitude and width of the wavepacket before and after interaction with a hydrodynamic state. It follows from the wave conservation law that the widths L ± of the wavepackets on both sides of the hydrodynamic state satisfy L − /k − = L + /k + (see (41)) so that Eqs. (44) and (45) yield the approximate conservation of amplitude: which agrees with Eq. (37) obtained for the limiting case of PW-mean flow interaction.

Wavepacket-rarefaction wave interaction
In this section, we consider in detail the interaction between a wavepacket and a RW; as we already mentioned, the linear wavepacket interacts with the RW only if it is initially placed at x = X 0 > 0 (see Tab. 1). The fields u and k are the solution of the Riemann problem studied in Sec. 3.2: the variation of u(x, t) is described by the relation (31) with V (u) = u, and the variation of k(x, t) is given in (39) obtained through the conservation of the adiabatic invariant (26). Although the linear wave is now localised, as we already mentioned, the definition of a wavenumber field k(x, t) for all x ∈ R still makes sense, the wavenumber of the carrier wave corresponds to k(x, t) evaluated at x = X(t), the position of the wavepacket. Since the wavepacket propagates with the group velocity v g (k, u), its position X(t) satisfies the characteristic equation dX dt = v g (k(X, t), u(X, t)) , X(0) = X 0 , with v g (k, u) = u − 3k 2 and some X 0 > 0. The integration of (47) is straightforward and yields where t + = X 0 /3k 2 + and t − = X 0 /3k + k − . Hence, during the interaction with the RW, the temporal variation of the dominant wavepacket wavenumber along the group velocity line is given by The wavepacket trajectory described by Eq. (48) is compared with the trajectory determined numerically in Fig. 6 for two different configurations: transmission (k + > k c , cf. (38)) and trapping (k + < k c ); the snapshots of the envelope a(x, t) of the wavepacket field ϕ(x, t) and the absolute value of its Fourier transform |φ(k, t)| are presented in Fig. 7. The numerical procedure implemented to extract ϕ(x, t) from the full solution for u(x, t) is explained in Appendix A. While the wavepacket propagates at constant velocity over a nonmodulated mean flow u(x, t) = u + or u(x, t) = u − , it gets decelerated during the propagation inside the RW for t + < t < t − ; note that acceleration here is understood as the increasing of the group speed |v g (X, t)| during  Figure 7: Numerical evolution of the wavepacket with an initial Gaussian envelope shape in the transmission interaction with a RW (left plot) and in the trapping interaction (right plot) (cf. Fig. 6). In both plots, the first column displays the extracted amplitude (envelope) a(x, t) of the wavepacket; the positions of the RW edges are shown by dotted lines ( ); the second column displays the amplitude of the Fourier transform of the wavepacket ϕ, denotedφ, which is initially also Gaussian. In the left plot, the wavepacket shape in Fourier space slightly deviates from Gaussian when it enters the leading RW edge (see the plots at t = 500); the wavepacket recovers its Gaussian form when it is fully inside the RW (t = 600) and after exiting the RW (t = 1300). Right plot: the wavepacket also deviates from its initial Gaussian form close to the leading edge of the RW (t = 500). Once inside the RW, the wavepacket recovers its Gaussian form and slowly decays as both its wavenumber and its amplitude go to 0.
propagation. If the transmission condition (38) is not satisfied, the incident wavepacket gets trapped inside the RW as its velocity converges asymptotically to the local background velocity (see Sec. 3.2). Moreover, as the wavelength of the carrier wave becomes infinitely large (K → 0), the amplitude of the wavepacket becomes infinitely small following the conservation of the integral of wave action (43), and the wavepacket eventually gets absorbed by the RW (see Fig. 7, right plot). We now draw certain parallels between the trapping of linear waves in RWs and the effect of "blocking" water waves in counter-propagating, inhomogeneous steady currents U (x), where the wavenumber k(x) also varies following the inhomogeneities of the current (recall the discussion in the Introduction). In this case, the adiabatic variation of the wavenumber is simply described by the conservation of the frequency ω(k(x); U (x)) = const and the wavenumber increases due to mean flow interaction. For waves near the critical, wave blocking wavenumber, the amplitude increases. In contrast to wave trapping due to wavepacket-RW interaction considered here, the wave blocking in the counter-propagating current is accompanied by a decrease in the wavelength of the wavepacket, until it reaches the stopping velocity v g = 0 at some finite wavenumber.
The trajectory of the wavepacket displayed in Fig. 6 shows that the wavepacket undergoes "refraction" due to its interaction with the RW. In the transmission configuration this results in both a speed shift and a "phase shift" of the transmitted wavepacket. The phase of the wave after its transmission is equal to . This result can also be obtained from the second adiabatic invariant pA in Eq. (33) where A = δk as in Eq. (40). Viewing the wavepacket as part of a fictitious periodic train of wavepackets, we recognise that the relative position of the wavepackets post-and-pre interaction X − /X 0 is inverse to the relative beating wavenumber shift Since k − < k + , the phase shift is negative in the considered situation. The formula (50) coincides with the result (41) obtained in the context of the wavepacket train. Fig. 8 displays the phase shift computed numerically for different wavenumber k + , which is in excellent agreement with relation (50).

Wavepacket-DSW interaction
We now consider the more complex case of linear wavepacket-DSW interaction. Such an interaction is generally described by two-phase KdV modulation theory, which is quite technical, with modulation equations given in terms of hyperelliptic integrals (Flaschka et al., 1980). The mean field approach adopted here enables us to circumvent these technicalities by employing the approximate modulation system (23) that yields simple and transparent analytic results that, as we will show, agree extremely well with direct numerical simulations. More broadly, the notion of hydrodynamic reciprocity described in Sec. 3.3 can be utilised without approximation to make specific predictions for wavepacket-DSW interaction for t > 0 based on wavepacket-RW interaction for t < 0. As already mentioned, wavepacket-DSW interaction admits two basic configurations (see Table  1): the transmission configuration, arising when X 0 > 0 and applicable to any incident wavenumber k + > 0, and the trapping configuration, when X 0 < 0 and the incident wavenumber k − > 0 is sufficiently small. The variation of u(x, t) inside the DSW is given by u = V −1 (x/t) (see Eq. (31)) where the characteristic velocity V (u) is defined by (7).  Fig. 6) and wavepacket-DSW interaction (right plot, cf., the trajectory in Fig. 11). In the two plots, the markers (pluses for the RW interaction and circles for the DSW interaction) correspond to the numerical simulation and the solid line to the analytical prediction (50). In full agreement with the theory, the phase shift is negative in the RW interaction (u − < u + ), whereas it is positive in the DSW interaction (u − > u + ).
The variation of the dominant wavenumber field k(x, t) in the wavepacket are given by the conservation of the Riemann (adiabatic) invariant q, which can be found by integrating the differential form Ξ (25). This differential form vanishes on the group velocity characteristic dx/dt = v g , yielding a relation between k and u specified by the ODE: with the boundary condition k(u + ) = k + if X 0 > 0 (k(u − ) = k − if X 0 < 0). Note that equation (51) for V (u) = u arises in the DSW fitting method where it determines the locus of the KdV DSW harmonic edge, see El (2005); . Here it has a different meaning and does not appear to be amenable to analytical solution because of the presence of elliptic integrals in the function V (u) (see Eq. (7)); in practice we solve (51) numerically. Once the relation k(u) has been determined, the (x, t)-dependence of the wavenumber inside the DSW is simply given by k = k(u(x, t)) = k(V −1 (x/t)). Finally, the trajectory of the wavepacket in (x, t) -plane is obtained by solving the equation of motion (47) with the fields u(x, t) and k(x, t) already determined. The results of our semi-analytical computations are presented in Figs. 9 and 10. Fig. 9 (left plot) displays the wave curves k(u) obtained from the numerical integration of (51) with (u − , u + ) = (1, 0). These wave curves can be seen as wavepacket trajectories in the parameter space (u, k) with the wavepacket dominant wavenumber K(t) = k(u(X/t)), where X(t) corresponds to the wavepacket position in physical space, and t being a parameter along the wave curve. The "propagation" of the wavepacket along the wave curve is then given by an ODE obtained by using u = V −1 (X/t), Eqs. (6), (47) and (51). Since the characteristic speed V (u) (6) of the Gurevich-Pitaevskii modulation equation is a decreasing function of u (cf. Fig. 2 shows that the dominant wavenumber of the wavepacket is increasing during its propagation inside DSW, in contrast to wavepacket-RW interaction for which V (u) = 1 > 0.
We now verify that the obtained approximate (mean-field) integral curves for the wavepacket-DSW interation are consistent with the transmission relation (36) for the linear PW-RW interaction, as expected from hydrodynamic reciprocity. In the transmission configuration, where k(u − ) > k c = 2 3 (u − − u + ) (see (38)), the wave curve k(u) is represented by a solid curve in Fig. 9 that connects u = u − to u = u + and, for a given incident wavenumber k + , the transmitted wavenumber is obtained by evaluating k(u − ). Fig. 9 right plot shows the comparison of the transmitted wavenumber k(u − ) evaluated by the above semi-analytical procedure with the value k − = k 2 + + 2 3 (u − − u + ) obtained from the wavepacket-RW transmission condition Eq. (36) by invoking hydrodynamic reciprocity. One can see the excellent agreement, which confirms the validity of the mean field approximation and its consistency with hydrodynamic reciprocity.
As expected, the behaviour of wave curves k(u) is drastically different for the trapping configuration, when k(u − ) < k c . In this case, the curves k(u) (represented by dashed curves in Fig. 9) do not connect u = u − to u = u + anymore, implying that the wavepacket initially placed at x = X 0 < 0 cannot reach the mean flow u = u + (trapping). Interestingly, these "trapping" wave curves k(u) are multi-valued, which implies that wavepackets with initial parameters u = u − , k < k c will return, asymptotically as t → ∞, to the DSW harmonic edge where u = u − . More specifically, the point u = u − , k = k c plays the role of an attractor in the parameter space for trapping con-  Fig. 9) with an initial phase-shift ∆(u + ) = 0; the dash-dotted curve ( ) corresponds to the limiting case k(u − ) 2/3. Right plot: deviation in the predicted phase shift relation ∆(u − ) from the actual phase shift ∆ − obtained from hydrodynamic reciprocity and eq. (50) as a function of ∆ − with u − = 1, u + = 0. |∆(u − ) − ∆ − |/∆ − is also plotted as a function of k − − k c , which can be mapped to ∆ − /X 0 through equation (50). same value k c = 2/3 with time. This "filtering" behaviour is unusual and drastically different from wavepacket-RW trapping (see Sec. 4.3), where the wavepacket trajectory is single-valued and k → 0 as t → ∞.
with q 0 = Q(k + , u + ). The last equality comes from the conservation of the adiabatic invariant P (q 0 , u)δk along the characteristicẊ = v g , cf. Sec. 2.2. Notice that, when the transmitted wavenumber k(u − ) is close to the critical value k c , the phase shift ∆(u − ) diverges (as illustrated by the dash-dotted line in Fig. 10). Similar to the transmission condition k − (k + ), the relation for the phase shift between the incident and the transmitted wavepacket can be recovered by evaluating ∆(u) at u = u − . The comparison between ∆(u − ) and the phase shift (50) is displayed in Fig 10  and shows good agreement between the two descriptions. Note that, similar to the determination of the transmission condition, the discrepancy |∆(u − ) − ∆ − | increases when k − gets closer to k c .
We now compare the wavepacket dynamics obtained through our modulation analysis (the resolution of (47), (51)) with the numerical solution of the KdV equation with initial conditions given by the partial Riemann data (42) (see Appendix A for details of the numerical procedure employed to trace the dynamics of a wavepacket inside a DSW). Fig. 11 displays trajectories for the two different configurations (transmitted wavepacket and trapped wavepacket), and snapshots of the envelope a(x, t) and the Fourier transform of ϕ(x, t) for the corresponding numerical simulation are presented in Fig. 12. The excellent agreement of the numerical simulations with the analytical predictions seen in Fig. 11 represents a further confirmation of mean field approximation employed in the derivation of the basic ODE (51).
Similar to the interaction with a RW, the group velocity of the linear wavepacket is not constant inside the DSW-but now the wavepacket accelerates in the transmission case and simultaneously experiences a wavenumber increase. Here, however, the variation of the dominant wavenumber K(t) along the trajectory is not everywhere adiabatic, see Fig. 11. This happens due to the violation of the basic inequality u x /u ϕ x /ϕ, which is at the heart of the modulation theory developed in Sec. 2.1 and used throughout this paper. The violation occurs due to the well-known logarithmic singularity of the mean flow gradient u x at the DSW soliton edge (Gurevich & Pitaevskii, 1974;El, 2005), see Fig. 2. The non-adiabatic behaviour of the wavenumber becomes obvious when one considers the evolution of the amplitude |φ(k, t)| of the Fourier transform of the linear field ϕ(x, t) along with the envelope a(x, t) of the field itself, see Fig. 12. Initially, in our simulations, both distributions have a Gaussian shape, but the Fourier transform amplitude distribution loses its unimodality when the wavepacket initially interacts with the leading, soliton, edge of the DSW (see Fig. 12, left plot at t = 1000). As a result, we are no longer in a position to define a nearly monochromatic carrier wave. Remarkably, however, as the linear wave progresses inside the DSW  Fig. 11), in the transmitted configuration (left plot) and in the trapped configuration (right plot). In both configurations, the first column displays the envelope amplitude of the wave-packet a(x, t) and the positions of the RW edges are indicated by dotted lines ( ); the second column displays the amplitude of the Fourier transform of the wavepacket ϕ, denotedφ. Left plot: close to the DSW leading edge, the wavepacket is composed of two different carrier waves with two different wavenumbers. Surprisingly, the wavepacket recovers its Gaussian form in physical space. Right plot: the determination of the wavepacket envelope is not straightforward and needs special care, the details are presented in appendix A.
the quasi-monochromatic structure of the wavepacket is restored (along with the Gaussian shape of the envelope).
The above described logarithmic divergence of the mean field gradient is absent in wavepacket-RW interactions considered in the previous section, where u = x/t inside the hydrodynamic state (see Eq. (31) with V (u) = u) such that u x /u remains finite but exhibits a discontinuity. We still observe a similar slightly non-adiabatic behaviour of the wavenumber in wavepacket-RW interaction at the initial stage, with u x /u not sufficiently small. This behaviour, expectedly, does not appear in the wavepacket -DSW trapping interaction (see Fig. 12, right panel), where the wavepacket, coming from the left of the hydrodynamic state, only interacts with the slowly varying part of the mean flow and there is no discontinuity in u x /u. Remarkably, Fig. 12 shows that the quasimonochromatic wavepacket structure is recovered when the interaction with the DSW edge is over so that its wavenumber is still described by the adiabatic analytical description when the wavepacket again propagates in the region where u x /u is small, and ultimately K = k − (where k − is given by the relation (36)) when the wavepacket is no longer interacting with the DSW. Note that this value of the wavenumber outside of the hydrodynamic state is expected to be in line with hydrodynamic reciprocity, since, as it was pointed out in Sec. 4.1, k is still a solution of a Riemann problem.
Finally, similar to wavepacket-RW interaction, we consider the phase shift of the wavepacket transmitted through the DSW. The numerical results are presented in Fig. 8 (right plot, lower panel) and are in very good agreement with the value (50) predicted analytically. One can also notice that for both RW and DSW interactions, the phase shift ∆ is given by the same relation (50), in accordance with the hydrodynamic reciprocity.
We note in conclusion that the trapping configuration is somewhat more difficult to treat analytically using the mean field approach employed in other cases. Although the trajectory, as well as the dominant wavenumber of the wavepacket, can be approximately described by our theory for short time evolution (see Fig. 11), we numerically observe that the dynamics of the DSW are no longer decoupled from the dynamics of the linear wave, and the distinction between the two structures becomes less and less pronounced after a sufficiently long time. The study of this problem is beyond the scope of the present work. We direct, however, the interested reader to Appendix A where the numerical determination of the linear wave for the trapping case is made for short time.

Conclusions and Outlook
We have introduced a general mathematical framework in which to study the interaction of linear wavepackets with unsteady nonlinear dispersive hydrodynamic states: rarefaction waves (RWs) and dispersive shock waves (DSWs). We use a combination of classical Whitham modulation theory and the mean field approximation to derive a new, extended modulation system that describes the dispersive dynamics of a linear wavepacket coupled to the nonlinear, long-wave dynamics of the mean flow in the hydrodynamic state. The mean field equation coincides with the long wave limit of the original dispersive equation when the hydrodynamic state is slowly varying (RW) but has a more complicated structure for rapidly oscillating states (DSWs) that requires a separate modulation analysis of nonlinear periodic waves. We show that the extended modulation system admits a convenient general diagonalisation procedure that reveals conserved adiabatic invariants during the course of wavepacket interaction with a slowly evolving mean flow. These adiabatic invariants predict transmission relations and trapping conditions for wavepackets interacting with hydrodynamic states. They also imply the general hydrodynamic reciprocity property whereby wavepacket interactions with RWs and DSWs exhibit the same transmission/trapping conditions. This enables one to circumvent a complicated analysis of the mean field behaviour in a DSW and to take advantage of the available wavepacket-RW relations to describe the transmission through a DSW or predict wavepacket trapping inside a DSW. This study has been performed using the KdV equation as a prototypical example, although the integrability properties of the KdV equation were not invoked, and the developed theory can be extended to other models supporting multi-scale nonlinear dispersive wave propagation.
The developed theory admits broad generalisations and opens a number of interesting perspectives. It can be naturally extended to physically relevant, non-integrable systems of "KdV type", such as the asymptotically equivalent, long wave Benjamin-Bona-Mahony equation, the Kawahara equation for capillary-gravity waves or the viscous fluid conduit equation. A number of immediate, intriguing questions arising in relation with the above systems are connections to possible effects of a non-convex linear dispersion relation, which is known to lead to profound effects on dispersive hydrodynamics, resulting in such novel phenomena as expansion shocks , DSW implosion (Lowman & Hoefer, 2013) and the existence of travelling DSWs (Sprenger & Hoefer, 2017). Another natural extension of this work is to the study of linear wave-mean flow interaction in the framework of integrable and non-integrable bidirectional systems such as the defocusing nonlinear Schrödinger equation or the Serre system for shallow water waves (Serre, 1953).
The "abstract" basic modulation system (1) has been extensively used in the theory of phase modulations revealing dispersive deformations arising near coalescing characteristics (see Bridges (2017); Ratliff & Bridges (2016) and references therein). At present, this theory does not include variations of the mean flow. The proposed here extension of the basic modulation system that couples modulations of the wavepacket with the mean field variations could also be useful for further development of the phase modulation theory.
The developed modulation theory of linear wave-mean flow interactions could be applied to the description of the interaction of wind-generated short waves with shallow water undular bores in coastal ocean environments as well as to the interactions of linear wavepackets with undular bores in internal waves. It can also be utilised in many physical contexts beyond classical fluid mechanics. In particular, similar to the soliton-mean flow interaction theory very recently developed in Maiden et al. (2018), it can be applied to a broad range of dispersive hydrodynamic systems describing wave propagation in nonlinear optics and condensed matter physics opening perspectives for experimental observation of the various interaction scenarios studied here. In fact, the linear wavepacket transmission and trapping configurations can be interpreted as hydrodynamic wavepacket scattering, a dispersive wave counterpart of hydrodynamic soliton tunnelling Maiden et al. (2018); Sprenger et al. (2018). In both cases, the role of a "barrier" or a "scatterer" is played by a large-scale evolving hydrodynamic state satisfying the same equation as the soliton (wavepacket). Finally, we mention the actively developing field of analogue gravity (see Barceló et al. (2011) and references therein), where the effects of dispersive wave trapping studied here may find interesting interpretations.  can be put in the form (3), we obtain the variations of the wavepacket ϕ(x, t) by evaluating the difference: ϕ(x, t) = u(x, t) − u H.S. (x, t) .
We then extract from the wavepacket amplitude a(x, t), the position of the wavepacket and from the spatial Fourier transformφ(k, t) = F[ϕ(x, t)], the wavepacket dominant wavenumber Note that here the position (57) and the dominant wavenumber (58) corresponds to average quantities instead of the maxima of a(k, t) andφ(k, t), respectively, which are not uniquely defined in some situations (cf. e.g. Fig. 12); when the wavepacket is Gaussian, the quantities (57), (58) are equivalent to the conventional definitions of the wavepacket position and dominant wavenumber.
The ansatz (3) proves to be inadequate to describe the variations of u(x, t) in the trapping interaction with a DSW (wavepacket coming from the left of the DSW, cf. Table 1). In fact, we observe in this case that the field u(x, t)−u H.S. (x, t) no longer corresponds to a quasi-monochromatic wavepacket, and exhibits, additionally to a Gaussian wavepacket, small harmonic excitations, cf. Fig. 13. We identify this deviation from the unimodality as a local phase shift of the DSW. Indeed, it is known that a soliton interacting with a dispersive wavetrain is phase shifted with respect to its free propagation (Ablowitz & Kodama, 1982), and a similar phenomenon could happen for DSWs which can be seen approximately as rank-ordered trains of solitons. Thus a schematic solution of the Riemann problem should read: u(x, t) = u H.S. (x − δ(x, t), t) + ϕ(x, t) , where δ(x, t) 1 corresponds to the small phase-shift induced by the wavepacket-DSW interaction. Yet we determine ϕ(x, t) numerically by computing the difference u(x, t) − u H.S. (x, t) such that: The shape of the oscillations of this new linear structure seems to correspond qualitatively to ∂ x u H.S. (x, t), cf. Fig. 13. The exact variation of the phase shift δ(x, t) is beyond the scope of the present work. Still, the evolution of the wavepacket can be recovered at early time (t ∼ 10 3 ), by eliminating in the Fourier transform of (56) the contribution clearly associated to the phase shift of the DSW. The resulting evolution is presented in Sec. 4.4, and comparisons to analytical predictions are displayed in Figs. 11 and 12. It is remarkable that, even if the dynamics of the DSW is slightly perturbed by the propagation of the wavepacket, the dynamics of the wavepacket is well described by the theory developed so far which confirms, once again, the mean field assumption for the wavepacket-DSW interaction.