Evolution of truncated and bent gravity wave solitons: the Mach expansion problem
Samuel Ryskamp, Michelle D. Maiden, Gino Biondini, Mark A. Hoefer
EEvolution of truncated and bent gravity wavesolitons: the Mach expansion problem
Samuel Ryskamp , Michelle D. Maiden , Gino Biondini , Mark A. Hoefer ∗ Department of Applied Mathematics, University of Colorado, Boulder, CO80309, USA, Department of Mathematics, State University of New York, Buffalo, NY, USA
Abstract
The dynamics of initially truncated and bent line solitons for theKadomtsev-Petviashvili (KPII) equation modelling internal and surfacegravity waves are analysed using modulation theory. In contrast to pre-vious studies on obliquely interacting solitons that develop from acuteincidence angles, this work focuses on initial value problems for the ob-tuse incidence of two or three partial line solitons, which propagate awayfrom one another. Despite counterpropagation, significant residual soli-ton interactions are observed with novel physical consequences. The initialvalue problem for a truncated line soliton—describing the emergence ofa quasi-one-dimensional soliton from a wide channel—is shown to be re-lated to the weak interaction of oblique solitons. Analytical descriptionsfor the development of weak and strong interactions are obtained in termsof interacting simple wave solutions of modulation equations for the localsoliton amplitude and slope. In the weak interaction case, the long-timeevolution of truncated and large obtuse angle solitons exhibits a decay-ing, parabolic wave profile with temporally increasing focal length thatasymptotes to a cylindrical Korteweg-de Vries soliton. In contrast, theresonant case of slightly obtuse interacting solitons evolves into a steady,one-dimensional line soliton with amplitude reduced by an amount pro-portional to the incidence slope. This strong interaction is identified withthe “Mach expansion” of a soliton with an expansive corner, contrastingwith the well-known Mach reflection of a soliton with a compressive cor-ner. Interestingly, the critical angles for Mach expansion and reflectionare the same. Numerical simulations of the KPII equation quantitativelysupport the analytical findings.
The oblique interaction of solitary waves or solitons is a fundamental problemin fluid dynamics and nonlinear sciences more broadly. Early theoretical con- ∗ Email address for correspondence: [email protected] a r X i v : . [ n li n . PS ] J u l igure 1: Contour plot of a line soliton solution and its slope parametrisation.The soliton propagation slope is q = tan ϕ .sideration of this problem for acute, collisional angles of incidence dates backto [32, 33] where weak and strong gravity water wave soliton interactions wereshown to be dependent upon the incidence angle and soliton amplitudes. In thecase of weakly interacting oblique solitons, a sufficiently small incidence angleleads to the approximate linear superposition of the two solitons accompaniedby a phase shift. The strong interaction case at large acute angles leads to a res-onant triad of obliquely interacting solitary waves, the so-called Miles resonantsoliton. Miles used his theory to identify the long-time dynamics of regular andMach reflection of a soliton incident upon a compressive corner or wedge.Building upon original water wave experiments by [37], Miles’ studies havesince been expanded and refined with the help of laboratory experiment [31,28, 24] and field observations [41, 1], numerical simulation [14, 39, 40, 9, 23, 28,19, 12] and exact N -soliton solutions [25, 7, 9, 11, 22, 24] of the Kadomtsev-Petviashvili equation [18] and its higher order generalisations( u t + uu x + u xxx ) x + u yy = 0 , ( x, y ) ∈ R , t > , (1)originally derived in the context of shallow water waves by [5] and for internalwaves by [16]. This version of the Kadomtsev-Petviashvili equation is known asthe KPII equation—the KPI equation occurs when + u yy → − u yy . The KPIIequation is a completely integrable equation [3] that admits a two-parameterfamily of stable line soliton solutions u ( x, y, t ) = a sech (cid:18)(cid:114) a
12 ( x + qy − ct ) (cid:19) , c = a q , (2)uniquely determined by the amplitude a > y -axis or slope q ∈ R . See Fig. 1 for a representative example.Although commonly identified with the oblique interaction of two shock waves[30, 26] or solitons [33], regular and Mach reflection can also be formulatedas an initial-boundary value problem in which a semi-infinite shock or soliton2igure 2: Schematic of Mach (left) and regular (right) reflection of a solitonimpinging upon a compressive corner.propagating parallel to a wall impinges upon a wedge or corner. Depending uponthe corner angle, the ensuing dynamics lead to the spontaneous generation ofa reflected wave and, in the case of irregular or Mach reflection, the additionalgeneration of a Mach or stem wave in which a resonant triad of three wavesmeet and propagate away from the wall. See Fig. 2 for a schematic of the tworeflection types.One approach to describe Mach reflection of solitons is the use of exact solu-tions of the KPII equation [28]. A classification of 2-soliton solutions in [25, 7, 11]was used to identify two particular 2-soliton solutions whose parameters can bechosen to satisfy the requisite structure of regular and Mach reflected waves. Todescribe the soliton-corner initial, boundary value problem, a nonlinear methodof images is applied and hypothesised to locally describe the long-time dynam-ics. This results in the critical angle ϕ cr —corner inclination measured from thepositive x -axis—for the transition from regular, ϕ > ϕ cr , to Mach, 0 < ϕ < ϕ cr ,reflection of an incident soliton with amplitude a astan ϕ cr = √ a. (3)Numerical simulations of initially V-shaped waves were used to justify the long-time, locally 2-soliton solution hypothesis for the soliton-corner initial, boundaryvalue problem [23, 19]. Herein lies a subtle difference between oblique solitoninteraction—described by exact 2-soliton solutions—and a soliton incident upona corner, which involves transient dynamics that, after long enough times, arelocally described by 2-soliton solutions.The aforementioned regular/Mach reflection problem involves a compressive corner with angle ϕ measured counterclockwise from the positive x -axis. Inthis paper, we consider the problem of a soliton incident upon an expansive corner, opening in the opposite direction so that ϕ is measured clockwise fromthe positive x -axis. This problem is rather different from the regular/Machreflection problem in many respects but we find an interesting parallel. Thecritical corner angle that separates regular expansion and Mach expansion foran incident soliton with amplitude a occurs precisely at ϕ cr , the same criticalangle separating regular and Mach reflection in Equation (3). Regular expansionoccurs when ϕ > ϕ cr and leads to the development of a decaying parabolic wavethat connects the incident soliton to the wall. Resonant or Mach expansion when0 < ϕ < ϕ cr involves the development of a new soliton perpendicular to the wall3igure 3: Left: diagram of Mach expansion for a front encountering a reversewedge boundary with small ϕ (top) and large ϕ (bottom). Note the emergenceof a new, smaller amplitude line soliton along the boundary when ϕ is small.Right: bent soliton initial conditions are identical for the purposes of analysis.with reduced amplitude relative to the incident soliton. The development of thisresonant soliton is the expansion analogue of the Mach stem in the reflectioncase. See Fig. 3 for a schematic of the two cases.A new approach is required to describe regular and Mach expansion of soli-tons because the dynamics are not described by 2-soliton solutions. While thetransient dynamics of Mach reflection are subtle, crucial aspects of Mach ex-pansion are transient and are not steady in a local reference frame. By makinga slowly varying assumption, we approximately describe the full expansion dy-namics developing from initial value problems by appropriate modulation of theline soliton (2). Our KPII numerical simulations demonstrate that modulationtheory is effective at uncovering key, quantitative features of the nonlinear wavedynamics.In order to describe the soliton-corner expansion problem, we take a seeminglycircuitous route by considering three classes of initial value problems, depictedin Fig. 4, for the KPII equation (1). From left to right, we identify the data astruncated, bent-stem, and bent soliton initial conditions. The bent soliton casecorresponds to the appropriate reflection, via the nonlinear method of images[24], needed to describe regular/Mach expansion as depicted in Fig. 3. Bentinitial conditions are equivalent to a single soliton moving along a boundary,where the boundary suddenly angles away from the front.The reason for considering these three classes of initial data in turn is bothmathematical and physical. From a mathematical point of view, each initialcondition limits to the next so that their solution is informed by the previ-4us and they share some solution properties. Moreover, these are among thesimplest and more natural kinds of non-solitonic initial conditions for the KPequation one could consider. Although the data consists of modulated solitons,it does not contain exact soliton solutions. In contrast to our approach, theprimary analytical means by which all previous studies have interpreted relatednon-solitonic initial value problems is by approximating the evolution of ini-tial configurations with exact N -soliton solutions of the KP equation. See, forexample, [9, 19, 12].From a physical point of view, the truncated soliton initial data models theemergence of a quasi-one-dimensional soliton from transverse confinement suchas the internal ocean solitons generated at the front of a river plume [36, 42] anda surface wave soliton when a channel suddenly widens. The rapid decelerationor transcritical propagation of a ship in open, shallow water can similarly launcha bent-stem or bent soliton from the ship’s prow [29], dependent upon theprow shape. Finally, internal wave solitons are ubiquitous in the world’s oceans[17, 41] and topography significantly impacts their propagation and interaction[43]. The classes of initial data in Fig. 4 represent cases where the modulatedsolitons propagate away from one another. Nevertheless, their interaction issignificant.We study these initial value problems using modulation theory in which thelocal soliton amplitude a = a ( y, t ) and slope q = q ( y, t ) are allowed to varyslowly in space and time. There are several approaches to derive the effectivemodulation equations. In Appendix A, we provide a derivation using multiplescale perturbation theory of the equations a t + 2 qa y + 43 aq y = 0 , (4a) q t + 2 qq y + 13 a y = 0 . (4b)The dynamical equation (4a) for the amplitude results from an appropriateorthogonality condition and the slope equation (4b) results from a consistencycondition of the modulated phase. The modulation equations (4) were alsoderived from a variable coefficient KP equation in [27] and by an averagedLagrangian approach in [35, 15]. They are also a limiting case of the moregeneral KP-Whitham modulation equations for periodic waves [2] in the caseof x independent modulations of a line soliton [8]. These soliton modulationequations are equivalent to the equations modelling the isentropic flow of apolytropic gas with density ∝ a / , velocity ∝ q , and ratio of specific heats γ = 5 /
3. The linearisation of the modulation equations (4) for small | q | wasused in Kadomtsev and Petviashvili’s original paper to determine the stabilityof line soliton solutions to the KP equation [18].Despite variation in only one spatial dimension, the corresponding modu-lated line solitons exhibit non-trivial two-dimensional structure. In particular,once a modulation solution a ( y, t ), q ( y, t ) is obtained, the modulated soliton is5igure 4: Initial data corresponding to a truncated line soliton (left), a bent-stem soliton (centre), and a bent soliton (right).reconstructed by projection onto (2) according to u ( x, y, t ) ∼ a ( y, t ) sech (cid:32)(cid:114) a ( y, t )12 ξ (cid:33) ,ξ = x + (cid:90) y q ( y (cid:48) , t ) d y (cid:48) − (cid:90) t c (0 , t (cid:48) ) d t (cid:48) , (5)where the soliton speed satisfies c ( y, t ) = a ( y, t ) / q ( y, t ) = − ξ t ( x, y, t ) (cf. (2)and (62)).The class of initial data we consider here corresponds to expansive condi-tions, so that we are guaranteed global existence of modulation solutions. Wewill obtain explicit modulation solutions in the form of simple waves and theirinteractions that describe the evolution of the data shown in Fig. 4.Our analysis is supported by numerical simulations of the KPII equation (1)using a Fourier pseudospectral method adapted from [19] that allows for out-going line solitons at the top and bottom of the simulation domain [ − L x , L x ] × [ − L y , L y ] through use of a windowing function. We maintain the nonlocal con-straint (cid:82) L x − L x u yy d x = 0 to high accuracy by including localised “image” initialdata whose superposition with the test data of interest satisfy this constraint.Simulations were terminated before the image and test data interacted. Forfurther details, see Appendix B. In this section, we summarise the classical analysis of the hyperbolic system (4).The modulation equations admit the amplitude symmetry a (cid:48) = a/A, q (cid:48) = q/ √ A, y (cid:48) = √ Ay, t (cid:48) = t, A > , (6a)6he quasi-rotational symmetry a (cid:48) = a, q (cid:48) = q + Q, y (cid:48) = y − Qt, t (cid:48) = t, Q ∈ R , (6b)the hydrodynamic scaling symmetry y (cid:48) = αy, t (cid:48) = αt, α > , (6c)and the reflection symmetry y (cid:48) = − y, q (cid:48) = − q, (6d)all leaving (4) unchanged in primed coordinates. Namely, if a ( y, t ) and q ( y, t )solve (4), so do a (cid:48) ( y (cid:48) , t (cid:48) ) and q (cid:48) ( y (cid:48) , t (cid:48) ).It was shown in [8] that appropriate linear combinations of equations (4a)and (4b) result in the equivalent pair of equations in characteristic form ± (cid:2) a t + (2 q ± √ a ) a y (cid:3) + 2 √ a (cid:2) q t + (2 q ± √ a ) q y (cid:3) = 0 , (7)which reveal the characteristic velocities U = 2 q − √ a and V = 2 q + √ a .Integration of (7) along each characteristic direction demonstrates that r = q − √ a, s = q + √ a (8)are Riemann invariants for the modulation equations (4), which afford the di-agonalisation r t + U r y = 0 , s t + V s y = 0 , (9a) U = 23 (2 r + s ) , V = 23 ( r + 2 s ) , (9b)where the characteristic velocities U and V are now written in terms of theRiemann variables r and s .Simple wave solutions of (9) correspond to variation in only one characteristicdirection so that one of the Riemann invariants r or s is constant, i.e., either q ( y, t ) + (cid:112) a ( y, t ) or q ( y, t ) − (cid:112) a ( y, t ) is independent of y and t . Since the char-acteristic velocities are ordered U ≤ V , we identify simple waves with variationalong the slow characteristics d x/ d t = U as 1-waves ( s = const ) and thosewith variation along the fast characteristics d x/ d t = V as 2-waves ( r = const ).Across a simple wave, the non-constant Riemann invariant’s characteristic ve-locity is monotonically increasing because the strictly hyperbolic system (4) isgenuinely nonlinear so long as a (cid:54) = 0 [8].Expansive initial data a ( y, q ( y,
0) corresponds to the condition that bothcharacteristic velocities U and V evaluated on the initial data are monotonicallyincreasing functions of y . By virtue of the fact that ∂U/∂r = ∂V /∂s = > r and s are non-decreasing functions of y .Simple waves propagate into constant regions of the y - t plane, and are there-fore fundamental building blocks for Riemann problems that posit step initial7ata at the origin. The interaction of simple waves can most conveniently beinvestigated by use of the hodograph transformation [13] in which the role ofdependent and independent variables is swapped. Namely, we take t = t ( r, s ), y = y ( r, s ), yielding the following set of linear equations t rs + 2 s − r ( t r − t s ) = 0 , (10a) y s = U t s , y r = V t r , (10b)so long as the Jacobian J = r t s y − s t r y remains nonzero. Equation (10a) isan Euler-Poisson-Darboux (EPD) equation that is equivalent to the radial waveequation in five dimensions, which admits the general solution t ( r, s ) = A + (cid:18) F ( r )( s − r ) (cid:19) r + (cid:18) G ( s )( s − r ) (cid:19) s , (11)for an arbitrary constant A ∈ R and functions F ( r ), G ( s ). In order to determine A , F , and G , one must specify suitable initial and/or boundary conditions. By a partial or truncated line soliton, we mean initial data in which the solitonmodulation amplitude is zero on one or two semi-infinite intervals, respectively.This scenario where a = 0 corresponds to a vacuum state in which the solitonslope q is undefined. To resolve this ambiguity, we will require that the valueof q in the neighbourhood of the point where a becomes positive correspondsto a simple wave in which one of r or s is constant [38]. Then, the propagationspeed of the vacuum front will necessarily be lim a → U = lim a → V = 2 q . This problem was previously posed and solved in [35] as a model for two-dimensional soliton diffraction. Here, we show that the resulting simple wavesolution forms an important building block for other, more complex, truncatedand partial soliton interactions.Without loss of generality (recall the symmetries (6)), we consider the partial(half) line soliton initial data a ( y,
0) = (cid:40) y > y ≤ , q ( y,
0) = 0 , y < , (12)for the modulation equations (4). This Riemann problem is solved by a 1-wave corresponding to a centred, self-similar rarefaction wave (1-RW) in which s = q + √ a ≡ U = y/t [35] (cid:113) a p ( y, t ) = t < y (cid:0) − y t (cid:1) − t < y < t y < − t , q p ( y, t ) = 1 − (cid:113) a p ( y, t ) . (13)8igure 5: Top: numerical simulation of partial soliton evolution for t ∈ (0 , , a → .
03 and a → .
97, respectively.These fronts are well approximated by straight lines with slopes given by thepredicted characteristic speeds from the solution in Eq. (13). Below, we will usethe solution in Eq. (13) to analyse the evolution resulting from more complexinitial conditions.
We now consider the KPII equation with initial conditions for a truncated linesoliton (recall Fig. 4 left) of length (cid:96) > a ( y,
0) = (cid:40) | y | ≤ (cid:96)/ | y | > (cid:96)/ , q ( y,
0) = 0 , | y | < (cid:96)/ . (14)By the reflection symmetry (6d), q and a are respectively odd and even func-tions of y for each t ≥
0. As in the case of the partial soliton, the truncatedsoliton slope in the vacuum region where a = 0 is determined by a simple wavecondition. Namely, for short times, a non-centred 1-RW is generated from theupper truncation point at y = (cid:96)/
2. Similarly, a 2-RW is generated from thelower truncation point y = − (cid:96)/
2. By use of the reflection symmetry (6d), theinitial evolution of these two simple waves can be represented as even or oddextensions of a shifted partial soliton (13) a ( y, t ) = a p ( | y | − (cid:96)/ , t ) , q ( y, t ) = sgn( y ) q p ( | y | − (cid:96)/ , t ) , (15)for y ∈ R . However, this solution only holds prior to the interaction of thesimple waves, which limits its validity to 0 ≤ t ≤ (cid:96) . A characteristic diagramshowing the 1-RW and 2-RW solutions emanating from Y = y/(cid:96) = ± / t = (cid:96) , the two simple waves intersect at y = 0. In order to understandwhat happens for T = t/(cid:96) > , we utilise the hodograph transformation and thecorresponding equations (10). The boundary conditions for the EPD equation(10a) can be obtained by recognising that, at the boundaries of the simple waveinteraction region, either r or s is constant. When s = 1, as in the 1-RWpropagating down from y = (cid:96)/
2, we differentiate the simple wave equation y − (cid:96)/ t = U = 23 (2 r + 1) (16)with respect to r to obtain the relation y r = 23 t r (2 r + 1) + 43 t. (17)Using this expression to eliminate y r from Eq. (10b), we obtain t r − − r t = 0 , s = 1 , r ∈ [ − , . (18a)Likewise, for the other 2-RW propagating up from y = − (cid:96)/
2, we obtain t s + 21 + s t = 0 , r = − , s ∈ [ − , . (18b)Integrating (18) from the initial time of simple wave interaction t ( − ,
1) = 3 (cid:96)/ t ( r,
1) = 3 (cid:96) (1 − r ) , r ∈ [ − , , (19a) t ( − , s ) = 3 (cid:96) (1 + s ) , s ∈ ( − , , (19b)10igure 6: Simple waves and their interaction in the characteristic Y - T ( y(cid:96) - t(cid:96) )plane for the truncated soliton initial data (14). The solid lines are simplewave characteristics. The shaded region corresponds to interacting simple wavesbounded by the dashed curves from Eq. (22).for Eq. (10a). Applying these boundary conditions to the general solution (11)yields A = 0, F ( r ) = (cid:96) (2 − (1 − r ) ), G ( s ) = − (cid:96) (2 − (1 + s ) ) and the solution t ( r, s ) = 3 (cid:96) (1 − rs )( s − r ) , ( r, s ) ∈ [ − , × ( − , . (20a)We now solve for y ( r, s ), by integrating both of Eq. (10b) y ( r, s ) = (cid:96) ( r + s )( r + 4 rs + s − r − s ) , ( r, s ) ∈ [ − , , (20b)where we used y ( − ,
1) = 0 at the initiation of simple wave interaction. Theexpressions (20) implicitly determine the simple wave interaction for t ≥ (cid:96) .We observe that the quantities Y = y(cid:96) , T = t(cid:96) (21)are independent of the truncated soliton length (cid:96) , a manifestation of the hy-drodynamic symmetry (6c). We will henceforth report results in the scaledvariables Y and T .The boundary of the simple wave interaction region is determined by eval-uating the hodograph solution (20) at s = 1 for Y > | Y | = 12 + 2 T − √ T , T ≥ . (22)These are the dashed curves in Fig. 6.As shown in Fig. 6 for short times, two non-centred simple waves describedby (13), (15) emanate from the soliton truncation points at Y = ± . For longtimes, the interaction boundary (22) approaches | Y | ∼ T with the same slope11s the outermost edges of the simple waves | Y | = 2 T + 1. Note, however thatthe two characteristic curves never cross.Returning to the physical variables a and q using (8), the simple wave inter-action is described by the hodograph solution (20), which yields the expressions Y = q a / (3 + a − q ) , (23a) T = 38 a / (1 + a − q ) . (23b)Since q (0 , T ) ≡ (cid:112) a (0 , T ) = 18 T (cid:18) f ( T ) / + 1 + 12 f ( T ) − / (cid:19) , ≤ T, (24a) f ( T ) = 18 + 12 T + T (cid:112) T . (24b)Using (23), one could obtain explicit expressions for a = a ( Y, T ) and q = q ( Y, T )for general Y and T . However, we can draw several important conclusions fromthe asymptotics of the implicit solution (23). For T (cid:29) | Y | (cid:28) T / , theamplitude is approximately independent of Y and the slope is approximatelylinear in Y a ( Y, T ) ∼ (cid:18) T (cid:19) / + 3 / (cid:18) T (cid:19) / ,q ( Y, T ) ∼ Y T + Y · / (cid:18) T (cid:19) / . (25)The one-term expansion for a and q in (25) is a self-similar solution of themodulation equations (4) [27].Using the above modulation solution to reconstruct the approximate soliton(2) yields interesting predictions for the initial data u ( x, y,
0) = (cid:40) sech (cid:16)(cid:113) x (cid:17) | y | ≤ (cid:96) | y | > (cid:96) (26)to the KPII equation (1). The soliton phase ξ = ξ ( x, y, t ) in Eq. (2) can beapproximated for large t using (25) as ξ ( x, y, t ) ∼ x + y t − (cid:18) (cid:96) (cid:19) / t / (27)Because the soliton maximum occurs where ξ = 0, we conclude that, for longtimes, the truncated line soliton shape approaches a moving parabola openingin the negative x direction with increasing focal length tx + y t = c ( t ) t, c ( t ) = (cid:18) (cid:96) t (cid:19) / . (28)12hile the parabolic shape is independent of the initial truncated soliton length (cid:96) , the wave speed is proportional to (cid:96) / . Concurrently, the soliton amplitudedecays according to Eq. (25).The shape and amplitude of the modulated line soliton within the simple waveinteraction region has an interesting connection to the cylindrical KdV (cKdV)equation. As noted in [4], introducing the change of variables u ( x, y, t ) = f ( η, t ) , η = x + y t (29)results in the exact reduction of the KPII equation (1) to the cKdV equation f t + f f η + f ηηη + 12 t f = 0 . (30)This equation admits slowly decaying soliton solutions [34]. Approximate solitonsolutions for t (cid:29) f ( η, t ) ∼ A ( t ) sech (cid:32)(cid:114) A ( t )12 ( η − z ( t )) (cid:33) , ˙ z ( t ) = A ( t )3 . (31)In order to determine the slowly varying amplitude A ( t ), we appeal to theconserved momentum P = (cid:90) R tf ( η, t ) d η (32)for any square integrable solution of cKdV (30). Inserting the slowly varyingsoliton ansatz (31) into (32), we obtain P = t A ( t ) / √ . (33)Thus, given some momentum P , the slowly varying amplitude is A ( t ) = (cid:18) / P t (cid:19) / . (34)If we choose the initial momentum to be P = 9 (cid:96) , then this amplitude equationmatches the leading order truncated soliton amplitude a ( Y, T ) (25) for large T .Moreover, the approximate cKdV soliton (31) admits the phase η − z ( t ) = x + y t − (cid:18) (cid:96) (cid:19) / t / , (35)which also matches the truncated soliton phase in Eq. (27), i.e., the leadingorder soliton slope q ( Y, T ) in (25).Figure 7 top depicts a numerical simulation of the truncated soliton initialdata (26) with (cid:96) = 300 that has been smoothed so as to minimise Gibbs typeoscillations [10]. Curved waves emanate from the truncation edges as the central13igure 7: Top: numerical evolution of truncated soliton initial data (70) accord-ing to the KPII equation for t ∈ (0 , , (cid:96) = 300. Bottom: modulatedsoliton amplitude a and slope q at noted times extracted from the numericalsimulation (solid curves) and the modulation solution (15), (23) (dashed curves)with the slightly different, fitted initial soliton length (cid:96) = 280. The dash-dottedblue lines correspond to the long-time asymptotic predictions (25) evaluated at t = 500. 14ortion propagates forward. When the central prominence decays, the entirewave forms a curved shape with decaying amplitude and curvature as timeincreases. These qualitative features are reflected in the obtained modulationsolution for the soliton amplitude a ( y, t ) and slope q ( y, t ) in Equations (15),(23). In order to quantitatively compare the simulation to modulation theorypredictions, we extract the modulated soliton amplitude and slope from thesimulation via a ( y, t ) = max x ∈ R u ( x, y, t ) , q ( y, t ) = − (cid:18) arg max x ∈ R u ( x, y, t ) (cid:19) y . (36)For the numerical computation of q , we smooth arg max u prior to differentia-tion. Figure 7 bottom displays the numerical (solid) and modulation (dashed)solutions. In order to quantitatively track the numerical simulation, we used theslightly smaller truncation width (cid:96) = 280 for the modulation solution in order toaccount for the smoothing of the initial data as given in (70). Both the solitonamplitude and slope closely match the full PDE evolution described by Eq. (1),demonstrating that our modulation analysis captures both the qualitative andquantitative features of the solution. The long-time ( T (cid:29)
1) asymptotic pre-dictions in Eq. (25) (dash-dotted) for a parabolic, decaying cKdV soliton alsocompare favourably with the numerical and modulation solutions for | y | (cid:47) (cid:96) despite the modest scaled time T = t/(cid:96) ≈ . The modulation solution for the truncated soliton consisting of two counter-propagating and then interacting simple waves motivates a broader class ofinitial conditions where we relax the assumption of zero soliton amplitude for | y | > (cid:96)/
2. In this section, we explore this scenario with three distinct config-urations: a special partially bent soliton, two bent solitons joined via a largeramplitude stem with nonzero (cid:96) , and finally a bent soliton with the same am-plitude throughout in which (cid:96) → We first consider a single bend at y = 0 where a ( y,
0) = (cid:40) y ≤ a y > , q ( y,
0) = (cid:40) y ≤ q y > , (37)For generic choices of 0 < a < < q <
1, the initial data (37) give riseto two separated simple wave solutions of the modulation equations (4): a fast15-RW and a slow 1-RW separated by a constant region. As a natural extensionof the partial line soliton solution (13), we restrict the data (37) so that only asingle simple wave—the 1-RW—is generated, i.e., s = q + √ a is constant: q + √ a = 1 , < a < , < q < . (38)We call the corresponding initial data (37) subject to (38) a partially bentsoliton, which will be useful to describe the bent-stem soliton initial data in thenext subsection. The constancy of s ( q ( y, t ) + (cid:112) a ( y, t ) = 1) and U ( a, q ) = y/t result in the 1-RW solution (cid:113) a pb ( y, t ) = a U t < y (cid:0) − y t (cid:1) − t < y < U t y < − t ,q pb ( y, t ) = 1 − (cid:113) a pb ( y, t ) , U = 2 − √ a . (39)This solution is the same as that of the partial soliton (13) for y < U t and limitsto the partial soliton modulation as the angled soliton amplitude vanishes a →
0. The positive amplitude, outgoing soliton gives rise to the slower characteristicvelocity U <
2. The evolution of a partially bent soliton (37) with a = 0 . q = 0 . U and − in Eq. (39) are favourably compared with the numerical simulation inFig. 8 bottom by identifying the front positions where a → .
67 and a → . We now consider the initial condition a ( y,
0) = (cid:40) | y | ≤ (cid:96)/ a | y | > (cid:96)/ , q ( y,
0) = (cid:40) | y | ≤ (cid:96)/ y ) q | y | > (cid:96)/ , (40)which is the modulation initial condition for the KPII data depicted in Fig. 4middle. This configuration describes an initial truncated soliton of length (cid:96) that is extended with outgoing line solitons of amplitude a < ± q . The case a = 0 and q = 1 corresponds to the truncatedsoliton (14). As similarly noted for the partially bent soliton, generic choices of0 < a < < q < y = ± (cid:96)/
2. The fastest and slowest waves, however,will not interact with the other waves, propagating far away from the initialstem region. These non-interacting, propagating waves are of less interest so werestrict the initial data such that a single simple wave is generated at each bend,as in the partially bent soliton case. Consequently, we assume the same simple16igure 8: Top: numerical evolution of the partially bent soliton according to theKPII equation for t ∈ (0 , , a = 0 . q = 0 .
3. Bottom: comparisonof the characteristic speeds of the upper (left panel) and lower (right panel) edgesof the partially bent soliton rarefaction wave. The plot displays the predictedfront speeds from the modulation solution (39) as reference lines (dash-dotted,blue), the numerically extracted front positions from the numerical simulation(solid, black), and a least squares linear fit (dashed, red) whose slope determinesthe measured speeds. 17ave constraint in Eq. (38) corresponding to a non-centred 1-RW emanatingfrom y = (cid:96)/ y = − (cid:96)/
2. We call thecorresponding initial data (40) a bent-stem soliton.This initial value problem is nearly identical to the truncated soliton problem.In fact, their solutions are essentially the same apart from one subtle yet crucialdifference: the velocities of the outermost edges of the counterpropagating sim-ple waves are different. These differing velocities lead to different interactionfeatures.We now use the partially bent soliton simple wave (39) to construct the coun-terpropagating simple waves for the bent-stem soliton initial data (40) a ( y, t ) = a pb ( | y | − (cid:96)/ , t ) , q ( y, t ) = sgn( y ) q pb ( | y | − (cid:96)/ , t ) , (41)for y ∈ R prior to simple wave interaction 0 ≤ t ≤ (cid:96) .Compared to the truncated soliton, the Riemann invariants for the bent-stemsoliton simple waves take values on the smaller square( r, s ) ∈ [ − , r ] × [ − r ,
1] where r = 1 − √ a < . (42)Consequently, the hodograph solution for the simple wave interaction region isthe same as for the truncated soliton, namely Eq. (20). However, the solutionmust be considered on the restricted domain (42). Two space-time characteristicdiagrams of the modulation solution for different values of a are shown in Fig. 9.The interaction region is shaded grey. The bottom point of the interaction regioncorresponds to the initiation of simple wave interaction when ( Y, T ) = (0 , ).Note that the characteristic for the uppermost edge of the incoming 1-RW, Y = + U T , eventually intersects the edge of the interaction region (22)at ( Y ∗ , T ∗ ). Similarly, the reflected characteristic emanating from Y = − intersects the interaction region at ( − Y ∗ , T ∗ ). These intersection points aregiven by Y ∗ = 12 + 3 − √ a a , T ∗ = 34 a (43)and are shown in the characteristic diagrams of Fig. 9 for two different choices of a . When a → T ∗ → ∞ and we recover the result for the truncated solitonin which the colliding simple waves do not completely intersect one another.For the bent-stem soliton in which 0 < a <
1, the existence of the intersectionpoints ( ± Y ∗ , T ∗ ) occurs because the characteristic velocity U is slower thanthe corresponding characteristic velocity of the truncated soliton U <
2. Thissubtle velocity difference leads to a significant change in the dynamics as wenow explain.For
T > T ∗ , the 2-RW that propagated from the lower bend at y = − (cid:96)/ r = r and expands along the upper, outgoing soliton. The uppermost, leading edgeportion of the simple wave is the straight-line characteristic Y = V ( T − T ∗ ) + Y ∗ , V = V ( r ,
1) = 23 (2 √ a + 1) . (44)18igure 9: Characteristic plots of interacting simple waves for the bent-stemsoliton initial data (40). Left: 0 < a ≤ , resulting in an infinite region ofinteraction for two simple waves in ( y, t )-plane. Right: < a <
1, resulting ina bounded interaction region. See main text for description.The boundary of the interaction region emanating from ( Y ∗ , T ∗ ) now becomesthe parametric curve Y = Y ( r , s ) , T = T ( r , s ) , s ∈ [ − r , , (45)where Y ∗ = Y ( r , T ∗ = T ( r ,
1) and the curve is traversed as s is decreasedfrom 1. A new Cauchy problem for the modulation equations (4) must be solvedwith data prescribed along the parametric curve (45). Because the region intowhich this Cauchy problem propagates is constant, ( a, q ) = ( a , q ), it is asimple wave, a 2-wave with r = r . The solution is determined by identifyingthe characteristics emanating from the boundary curve (45). Given any s ∈ (max( − r , ,
1) along the boundary curve (45), the corresponding characteristicalong which s is constant is the straight line Y = 23 ( r + 2 s )( T − T ( r , s )) + Y ( r , s ) . (46)Example characteristics are shown in Fig. 9 right.A bifurcation occurs in the shape of the interaction region depending on theinitial outgoing soliton amplitude a . For sufficiently large a , the interactionboundary (45) terminates when Y = 0, which from the hodograph solution(20b) occurs when s = − r . According to the parametric curve (45), it wouldappear that s = − r can occur for any − ≤ r <
1. However, the hodographsolution for T in Eq. (20a) shows that T → ∞ as s → r . As s is decreasedfrom 1 in the parametric curve (45), s attains the value r before it reaches − r if and only if r >
0. Consequently, the critical value r = 0 determinesthe bifurcation from an unbounded (when 0 ≤ r ≤
1) to a bounded (when − < r <
0) simple wave interaction region. We now consider each case inturn.The characteristic diagram for an unbounded interaction case where r < < a ≤ / / ≤ q <
1) is shown in19imple wave interaction region unbounded boundedconstraints on initial data 0 ≤ a ≤ ⇐⇒ ≤ q ≤ < a < ⇐⇒ < q < initial geometric constraints strongly bent soliton ( q > a ) weakly bent soliton ( q < a )long time dynamics decaying parabolic soliton non-decaying line solitonTable 1: Dynamics of the bent-stem soliton initial data (40).Fig. 9 left. Aside from the intersecting characteristics at ( ± Y ∗ , T ∗ ) and theconcomitant simple wave (46) that emerges from the interaction region, thebent-stem soliton diagram is similar to the truncated soliton solution in which a = 0 (cf. Fig. 6). In fact, the long-time asymptotic behaviour of the solutionis identical to the truncated line soliton (25) when | Y | ≤ T / in which the stemforms a decaying soliton that approaches the parabolic-shaped cKdV soliton(31).In contrast, when r < / < a < < q < / Y c = 0. Since s = − r determines theterminus of interaction, the corresponding closing time T c can be calculatedfrom the hodograph solution (20a) T c = − r r = 3(1 − √ a + 2 a )4(2 √ a − . (47)The corresponding soliton slope is zero by reflection symmetry and the ampli-tude can be read off from s = − r , giving a ( Y c , T c ) = a c = (2 √ a − , q ( Y c , T c ) = q c = 0 . (48)From this closing point, a constant region emerges, bounded by the edges of thesimple wave (46) and its symmetric reflection | Y | = V ( r , − r )( T − T c ) = − r ( T − T c ) = 23 (2 √ a − T − T c ) . (49)This constant region corresponds to the emergence of a line soliton with am-plitude 0 < a c <
1. Our findings for the bent-stem soliton are summarised inTable 1.Figure 10 depicts the numerical evolution of bent-stem solitons for each ofthe scenarios in Table 1. For the panels in Figure 10 top, the initial conditionsare nominally √ a = 0 . q = 0 .
2, with (cid:96) = 100. From our analysis, theemergence of a constant region in the modulation, i.e., a vertical soliton withamplitude a c = 0 .
36, should begin to appear at t = (cid:96)T c ≈ t = 400 inFig. 10 top, right, the vertical soliton has emerged with amplitude very closeto the predicted value 0 .
36 shown in Fig. 11 left. In contrast, for the panels in20igure 10: Top: numerical simulation of the bent-stem soliton for √ a = 0 . q = 0 . t ∈ (0 , , t = (cid:96)T c = 236 (cf. Fig. 9right). Bottom: numerical simulation of the bent-stem soliton for √ a = 0 . q = 0 . t ∈ (0 , , (cid:96) = 100. 21igure 11: Comparison between numerical simulation (solid line) and the mod-ulation solution (dashed line) of the bent-stem amplitude decay at y = 0 for theparameters in Fig. 10 top (left) and bottom (right). In order to account for thesmooth initial data, the modulation solution is “fitted” by choosing (cid:96) = 80.Figure 10 bottom, the initial conditions are √ a = 0 . q = 0 .
7, again with (cid:96) = 100. As expected, the system forms a parabola which slowly decays overtime.For quantitative analysis, we consider the amplitude decay at y = 0 for thebent stem simulations in Figure 11. On the left is displayed the amplitude decayfor the weakly bent simulation shown in Fig. 10 top, while on the right is datafrom the strongly bent simulation from Fig. 10 bottom. Note especially thelong-time behaviour. For the weakly bent stem, the amplitude asymptoticallyapproaches a c = 0 .
36 as predicted, while for the strongly bent stem case, theamplitude continues to decrease as t → ∞ . As in the truncated case, we slightlyreduce the length (cid:96) in the modulation solution to (cid:96) = 80 in order to account forthe smoothing of the initial conditions.We also consider the predicted soliton phase ξ compared to the numericalsimulations. This is shown in Figure 12. The overlaid predicted phases (dashedcurves) were generated by using the modulation solution for q and then numer-ically integrating for ξ according to Eq. (5). We utilised the speed c ( y, t ) in theprediction after fitting the phase so that it lines up with the front’s maximumalong y = 0 at t = 100 in the leftmost panels. The ensuing phase profiles at t = 200 and t = 400 are slightly advanced relative to the numerical simulation,which can be attributed to higher order phase errors that are common in solitonperturbation theory. Importantly, the shape of the front’s crest is well-describedby the modulation solution for q .It is evident that both the weakly and strongly bent-stem numerical evolutionsare well approximated by the modulation solution, asymptoting to a line solitonand a parabolic wave in long time, respectively.22igure 12: Modulation solution phase (dashed) overlaid on contour plots forweakly (top) and strongly (bottom) bent-stem initial conditions when t ∈ (100 , , .3 Bent solitons We now consider the bent-stem soliton initial data (40) with a vanishing stem (cid:96) →
0, i.e., a bent soliton in which a ( y,
0) = a , q ( y,
0) = (cid:40) q y > − q y ≤ . (50)Figure 4 right displays the corresponding initial condition for the KPII equation(1). In contrast to the truncated and bent-stem soliton initial conditions, theinitial conditions (50) for the modulation equations (4) correspond to a Riemannproblem. We limit our consideration to an expansive Riemann problem bytaking q >
0. This case corresponds to a partial soliton interacting with anexpansive corner (cf. Fig. 3). If q <
0, a case we do not consider, the partialsoliton interacts with a compressive corner and gives rise to regular and Machreflection (cf. Fig. 2).For the special case where q + √ a = 1, cf. (38), the bent soliton’s evolutioncan be obtained directly from the bent-stem soliton evolution by taking the van-ishing stem limit (cid:96) →
0. Consequently, the bent soliton inherits the bent-stemsoliton’s bifurcation in long-time dynamics. When √ a > q , the modulationsolution for the bent-stem soliton post simple wave interaction ( T > T c ) exhibitsan expanding constant region with a = a c , q = q c in Eq. (48) bounded by thecharacteristics (49). We refer to this as the strong interaction case.When √ a < q and q + √ a = 1, we can take the (cid:96) → ± y ∗ , t ∗ ) = ( ± (cid:96)Y ∗ , (cid:96)T ∗ ) → (0 , y = 0. However,the characteristics (46) leaving the boundaries of the interaction region persist.At y = 0, the soliton amplitude in the interaction region is explicitly (24a) sothat a (0 , t ) → (cid:96) →
0. This case corresponds to weak interaction.In order to elucidate more details, generalise the aforementioned results, andprovide an alternative method of solution, we now solve the Riemann problem(50) for the bent soliton directly. We relax the assumption q + √ a = 1 andconsider general a > q > √ a > q . The upper andlower solitons cannot be connected by a single simple wave, which would require √ a − q = √ a + q , leading to the conclusion q = 0. Instead, we introducethe intermediate state ( a i , q i ) and connect it to ( a , ± q ) with simple wavessatisfying √ a i − q i = √ a − q , √ a i + q i = √ a + ( − q ) . (51)The solution to these equations under the given constraints is q i = 0 , √ a i = √ a − q , (52)where we use a 2-RW to connect the intermediate state to the top soliton anda 1-RW to connect to the bottom soliton.24igure 13: Numerical simulation of bent solitons for the strong interaction when a = 1, q = 0 . t ∈ (0 , , a = 1, q = 1 . t ∈ (0 , , (cid:112) a ( y, t ) = √ a V ( a , q ) t < | y | ( yt + 2 √ a i ) V ( a i , t < | y | < V ( a , q ) t √ a i | y | < V ( a i , t ,q ( y, t ) = sgn( y ) q V ( a , q ) t < | y | (cid:112) a ( y, t ) − √ a i V ( a i , t < | y | < V ( a , q ) t | y | < V ( a i , t , (53)where V ( a, q ) = 2 q + √ a . The solution contains a vertical line soliton expand-ing in y with amplitude a i in Eq. (52). This solution agrees with our analysis ofthe (cid:96) → q + √ a = 1. As we will show inthe next subsection, this strong interaction case corresponds to Mach expansionof a soliton interacting with a corner.For weak interaction when √ a < q , the above calculation fails because √ a i q cr = √ a (54)between the two classes of bent soliton dynamics. Instead, we introduce theintermediate vacuum state a i = 0. In order to connect to vacuum with a simplewave from each of the bent solitons, we require q i (cid:54) = 0. Since q in the vacuumregion is undefined, we determine it locally based on the simple wave criterion.By symmetry, q i+ = lim y → + q ( y, t ) = − lim y → − q ( y, t ) = q i − , (55)with the initial discontinuity at y = 0. We can use Riemann invariants tocalculate the values of q i ± . For the top simple wave, √ a − q = − q i+ , and bysymmetry q i+ = − q i − . The solution for the top simple wave is (cid:112) a ( y, t ) = √ a V ( a , q ) t < | y | ( yt − q i+ ) V (0 , q i+ ) t < | y | < V ( a , q ) t | y | < V (0 , q i+ ) t ,q ( y, t ) = sgn( y ) q V ( a , q ) t < | y | (cid:112) a ( y, t ) + q i+ V (0 , q i+ ) t < | y | < V ( a , q ) tq i+ | y | < V (0 , q i+ ) t , (56)with a symmetric reflection for the bottom simple wave. This solution consists ofan expanding vacuum a = 0 region connected to outgoing, canted line solitonsby simple waves. These simple waves are rotated versions (see Eq. (6b)) ofthe partial soliton solution (13). They are completely disconnected from oneanother; at this order of approximation, the interaction is negligible in that theevolution of the upper and lower branches can be analysed independently of oneanother. The vacuum region’s rate of expansion is proportional to q −√ a > q > q cr ) causes the outgoing solitons to separatefrom one another sufficiently fast so that their interaction is negligible withinthe context of modulation theory. As we will show in the next subsection, thisweak interaction case corresponds to regular expansion of a soliton interactingwith a corner.The numerical simulations in Figure 13 are essentially consistent with thesepredictions. For the bent soliton with a = 1 and q = 1 . > q cr = 1 in Figure 13bottom, a decaying parabolic front with trailing oscillations appear. Althoughan expanding, strictly vacuum region predicted by modulation theory is notimmediately apparent, amplitude decay is present. We have verified that theamplitude, shape, and propagation of the leading parabolic front is consistentwith the profile for the cKdV parabolic soliton (31). The trailing oscillationsare consistent with a two-dimensional generalisation of the oscillatory shelf thatis common for perturbed KdV (cKdV) problems. In contrast, for an initial bentsoliton with a = 1 and q = 0 . < q cr = 1, as seen in Figure 13 top, a newvertical line soliton with reduced amplitude appears.26igure 14: Modulated soliton amplitude a and slope q extracted from the nu-merical simulation in the strong interaction case of Fig. 13 top (solid curves)and the modulation solution (53) (dashed curves) at different times.Quantitative results further confirm our analysis. In Figure 14, it is evidentthat the predicted solution for both soliton amplitude and slope, extracted ac-cording to Eq. (36) captures the behaviour of the strong interaction case with a = 1 and q = 0 . < q cr = 1. As expected, the solution for large timesapproaches a line soliton with a i ≈ . q = 1 . > q cr = 1 shown in Fig. 15 ( (cid:96) = 0 case),the simulation’s lead wave slope (solid) is well approximated by the modulationsolution (dash-dotted). The amplitude does not reach zero as modulation the-ory predicts, although it does continually decrease. Recalling that modulationtheory applies under slowly varying assumptions, it is not surprising that animmediate transition from unit amplitude to zero amplitude does not occur inthe numerical simulation. In Figures 11 and 14, we observe that the numericalsolution temporally lags behind the modulation solution. The same happenshere in Fig. 15, albeit to a more significant degree in amplitude. While thisweakly interacting bent soliton simulation deviates from the modulation solu-tion in amplitude, in fact it can be reasonably approximated by the bent-stemmodulation solution for moderate stem length (cid:96) . This is to be expected; thebent-stem analysis for sufficiently large q shows that any positive stem length (cid:96) > (cid:96) > (cid:96) = 12 to the numericalsimulation of the weakly interacting bent soliton in Fig. 15. Now the decayingparabolic front is represented in the modulation solution. In other words, theweakly interacting bent soliton evolution exhibits a remnant of the bent-stemsoliton solution. 27igure 15: Comparison of bent-stem soliton modulation with small stem (cid:96) = 12(dashed, red) and no stem (cid:96) → We are now in a position to interpret this analysis in the context of the soliton-corner initial boundary value problem that is schematically depicted in Fig. 3.Consider a vertical, partial soliton with amplitude a propagating in the positive x direction adjacent to a horizontal wall located at y = 0. When this solitonencounters a corner at the origin that suddenly opens or turns away by theclockwise angle ϕ >
0, it expands. The nature of its expansion depends on thecorner angle and soliton amplitude. Via the nonlinear method of images, we mapthe partial soliton at the moment it encounters the corner to the bent solitoninitial data for the numerical simulation in Fig. 4 right and for modulationtheory in Eq. (50) with a = a and q = tan ϕ . Consequently, the critical cornerangle separating two distinct types of soliton-corner interaction is (cf. Eq. (54))tan ϕ cr = √ a. (57)For the case ϕ > ϕ cr = arctan √ a (sharp corner), the soliton almost com-pletely separates from the wall. The residual soliton-wall interaction is througha decaying parabolic soliton. In turn, the propagating partial line soliton de-cays, retreating further away from the wall. We term this case regular expan-sion . In contrast, for slight turns of the wall at the corner where 0 < ϕ < ϕ cr =arctan √ a , the soliton also develops a curved front that instead terminates at anon-decaying soliton perpendicular to the wall with lower amplitude than theincident soliton. The predicted soliton wall amplitude is a w = ( √ a − tan ϕ ) . (58)Despite propagating away from the wall, it does not “escape” its influence likein the regular expansion case. The residual soliton formed at the wall is theexpansion counterpart to the Mach stem that forms during the course of Machreflection (cf. Fig. 2). We term this case Mach expansion . Surprisingly, the28oliton-corner expansion regular expansion Mach expansionsoliton amplitude a , corner angle ϕ < √ a < tan ϕ tan ϕ ≤ √ a long time dynamics at wall decaying parabolic soliton non-decaying line solitonTable 2: Soliton-expansive corner initial, boundary value problem.crossover from regular to Mach expansion occurs at precisely the same cornerangle (57) as the crossover from regular to Mach reflection (3). These resultsare highlighted in Table 2. Using the KPII equation as a model of multidimensional gravity wave solitons,we describe the evolution of truncated and interacting oblique solitons usingmodulation theory, and we compare the analytical predictions with the resultsfrom direct numerical simulations. The initial value problems considered aredistinguished by geometric configurations of partial solitons that propagate awayfrom one another. Despite this, residual interactions between solitons occur thatlead to nontrivial wave patterns.An initial soliton that is transversely confined or truncated, and that propa-gates into an open region, first “curls” at the endpoints as it then morphs intoa parabolic shape over long time. The front’s parabolic shape flattens with alinearly increasing focal length with time. The parabolic wave’s amplitude andspeed exhibits algebraic decay proportional to t − / . Such wave patterns ap-pear to be common in images of oceanic internal waves for near-shore conditions[17, 36, 42].We also generalise the above truncated soliton configuration by appendingcanted partial solitons to it and find new dynamical behaviour. In addition tothe decaying parabolic wavefront for sufficiently canted solitons, a non-decayingvertical soliton with reduced amplitude relative to the original soliton segmentappears. This bifurcation in behaviour carries over to the soliton-corner expan-sion problem.The final initial value problem for a bent soliton also describes the interactionof a soliton propagating parallel to a wall with an expansive corner. For a sharpenough bend, the solitons exhibit weak interaction through a decaying parabolicfront, the case we identify as regular expansion of a soliton. For a slight enoughbend, the solitons continue to interact in a resonant way so as to produce anintermediate, non-decaying soliton at the wall with reduced amplitude a w (58)that connects them. This case of Mach expansion parallels the well-knownMach reflection of oblique solitons, both of which occur at the same criticalangle and display a similar transition between strong and weak interactions.Such a transition in the wave dynamics could potentially be observed in theshallow water context by soliton generation from a moving disturbance [27, 29]29r by experiments analogous to previous shallow water studies involving Machreflection [37, 31, 28, 24].In the Mach reflection case, an important quantitative test of the theory isits prediction of amplitude amplification of the Mach stem at the wall. Afterproperly taking into account a higher-order asymptotic approximation to shal-low water waves than the KP equation, the soliton wall amplitude has beendemonstrated to satisfactorily predict experiments across a range of incidentsoliton parameters (angle and amplitude) [24]. This points to a possible quanti-tative test of the Mach expansion theory presented here by measuring the wallreduction factor α ≡ a w a = (cid:18) − tan ϕ √ a (cid:19) , tan ϕ < √ a (59)in similar shallow water experiments. The factor α < a w ) and pre ( a ) corner interaction, respectively.These results also motivate the conjecture that outgoing gravity line solitonspropagating away from one another with slopes ± q ∞ and similar amplitudes a ∞ lead to a decaying parabolic or negligible interaction region when sufficientlysloped q ∞ ≥ √ a ∞ but leave a residual line soliton between them when q ∞ < √ a ∞ .This work demonstrates the practical utility and efficacy of soliton modula-tion theory to describe rich nonlinear wave dynamics. All the solutions that weconsider are globally existing simple wave or interacting simple wave solutionsof the hyperbolic modulation equations. These solutions, when projected backonto a line soliton, quantitatively agree with direct numerical simulations ofthe KPII equation. Although not previously recognised as such, simple wave-modulated solitons can also be seen in a variety of previous KPII numericalstudies [14, 23, 19, 12]. In particular, the transient portion of the reflected wavethat develops during both regular and Mach reflection of a soliton by a cornerappears to show a similar wave pattern to the partial soliton studied here. Anintriguing problem is to consider the modulation equations with initial data thatis compressive, i.e., that would give rise to shock solutions. Indeed, two collid-ing partial solitons and “V-shaped” initial conditions [23, 12] for regular andMach reflection, give rise to compressive Riemann problems for the modulationequations (4). How are such initial value problems regularised? What do shocksmean in this modulation context? Acknowledgements
The work of MAH and SR was supported by NSF grant DMS-1816934. Thework of MM was supported by the NSF GRFP. The work of GB was supportedby NSF grant DMS-2009487. Authors thank the Fields Institute Focus Programon Nonlinear Dispersive Partial Differential Equations and Inverse Scattering inthe summer of 2017 where this research was initiated.30 eclaration of interests
The authors report no conflict of interest.
A Derivation of the soliton modulation equa-tions
Here we show how Eqs. (4) can be directly derived from the KP equation (1)using multiple scales, without employing the full Whitham modulation theory.First, we introduce the rescaling X = (cid:15)x, Y = (cid:15)y, T = (cid:15)t (60)into the KPII equation (1) (cid:0) u T + uu X + (cid:15) u XXX (cid:1) X + u Y Y = 0 . (61)In order to study modulated line solitons, the asymptotic expansion u ( X, Y, T ; (cid:15) ) = u ( ξ, Y, T ) + (cid:15)u ( ξ, Y, T ) + · · · , u ( ξ, Y, T ) = a ( Y, T )sech ( ηξ ) ,ηξ X = 1 (cid:15) (cid:114) a ( Y, T )12 , ηξ Y = 1 (cid:15) (cid:114) a ( Y, T )12 q, ηξ T = − (cid:15) (cid:114) a ( Y, T )12 c ( a, q ) , (62)is assumed where ξ is the fast variable. The coefficient η ( Y, T ) = (cid:112) a ( Y, T ) / ξ XY = ξ Y X = 0. The consistencycondition ξ Y T = ξ T Y yields the slope modulation equation (4b). The amplitudemodulation equation (4a) is obtained by inserting the ansatz (62) into the KPIIequation (61). At first order in (cid:15) , we obtain an inhomogeneous ODE for u that,when integrated once with respect to ξ , is − c∂ ξ u + ∂ ξ ( u u ) + ∂ ξξξ u + q ∂ ξ u = − ( ∂ T u + 2 q∂ Y u + q Y u ) . (63)Solvability over the space of L ( R ) solutions is enforced by the orthogonalitycondition( ∂ T + 2 q∂ Y + q Y ) (cid:90) R u d ξ = ( ∂ T + 2 q∂ Y + q Y ) (cid:32) √ a / (cid:33) = 0 , (64)which, upon simplification, results in the amplitude modulation equation (4a). B Numerical integration of the KP equation
To validate our analytical results, we implement the pseudospectral method de-scribed in [19], which utilises a hyper-Gaussian windowing function and Fourier31iscretisation for non-periodic data in y and periodic data in x . We essentiallyfollow [19] with a few modifications. The method proceeds as follows. Insteadof solving the KP equation (1) for the original function u , we instead solve anequivalent PDE for a windowed function v ( x, y, t ) = W ( y ) u ( x, y, t ) , (65)where W ( y ) = e − a n | y/L y | n , (66)with n = 27 and a n = 1 . n ln 10. The rapid decay of W ( y ) for | y | near L y ensures that v = 0 near the top and bottom boundaries of the domain so thatits periodic extension is smooth. In this region, we assume that the solution u asymptotes to non-modulated line solitons. Thus, u can be decomposed as u = v + (1 − W )ˆ u, (67)where ˆ u are line solitons with constant a and q of the form shown in Eq. (2).Inserting the transformation (67) into the KP equation (1), we obtain an equiv-alent PDE for the windowed function v ( v t + vv x + v xxx ) x + v yy = (1 − W )( W ˆ u ˆ u x − ( v ˆ u ) x + 2 W (cid:48) ˆ u y + W (cid:48)(cid:48) ˆ u ) x , (68a)subject to the initial conditions v ( x, y,
0) = u ( x, y, − (1 − W ( y ))ˆ u ( x, y, . (68b)The above PDE (68) is solved numerically. At each time t we reconstruct thetrue solution u using Eq. (67). The advantage of this method is that since v iszero at the domain boundaries, we obtain spectral convergence using a Fourierdiscretisation in space. In order to preserve spectral accuracy, the derivatives ˆ u x ,ˆ u y , W (cid:48) , and W (cid:48)(cid:48) on the right hand side of (68a) are calculated analytically; thisis one difference from [19] where these derivatives are calculated using finite dif-ference approximations. Time stepping is performed with an integrating factorand the classic fourth-order Runge-Kutta scheme. Simulations are terminatedbefore the windowing region is corrupted by non-solitonic data.The numerical scheme described above is validated using an exact Y-shapedsolution, also known as the Miles resonant soliton (see Fig. 16 left). By refiningthe time and space steps, convergence is obtained at t = 10 to approximately10 − in the 2-norm relative difference of the numerical and exact travelling wavesolutions, as shown in Figure 16. The rate of spatial convergence in Figure 16demonstrates that spectral accuracy is obtained. In order to ensure reasonablecomputation time and memory demands, the simulation parameters are fixedat the grid spacing ∆ x = 1 / t = 10 − . Based on Figure 16,this ensures 2-norm errors below 10 − for resolving the Y-soliton solution upto t = 100. Simulations are run on domains of various sizes that depend uponthe problem, typically with an area comparable to [ − , . To take fulladvantage of our system’s graphics processing unit, we used single precision forall calculations presented in the main text.32igure 16: Convergence of Y-shaped soliton (left) in space (centre) and time.These simulations were run on the domain [ − , × [ − ,
64] with ∆ x = ∆ y .For the centre plot, we fixed ∆ t = 10 − , and for the right we fixed ∆ x = ∆ y =1 / (cid:82) u yy dx = 0 (see [6, 20]). The bent solitoninitial conditions satisfy this constraint. For the initial conditions which donot satisfy the constraint, a reflection is added that ensures that (cid:82) u y dx = 0,thereby satisfying the constraint. The reflected solitons have parameters ( a r , q r )defined by a r ( y ) = (1 − (cid:112) a ( y, , q r = 0 . (69)The full initial conditions including the reflection for the truncated and bent-stem cases are displayed in Fig. 17. Simulations are terminated before thereflected solitons influence the solution in the region of interest.To reduce Gibbs phenomenon, initial data for the simulations are obtainedby smoothing the discontinuous initial data a ( y,
0) and q ( y,
0) and inserting thisdata into Eq. (5). For example, the truncated soliton data (14) becomes u ( x, y,
0) = a ( y, (cid:32)(cid:114) a ( y, x (cid:33) , a ( y,
0) = 12 (cid:18) tanh (cid:18) (cid:96)/ − | y | w (cid:19) + 1 (cid:19) , (70)with w = 5, (cid:96) = 300.The two parameters of interest in this paper are the soliton amplitude and thesoliton centre. The amplitude is obtained from the numerical results by findingthe maximum value over x ∈ [ − L x , L x ] for each fixed y and t in the domain vialocal interpolation of the solution. The soliton centre is simply the location ofthat maximum value. These are then compared to the predicted amplitude andcentre values. The predicted soliton centre is calculated using the analyticalsolution for q ( y, t ) combined with (5) or, if y is in the simple wave interactionregion, using (28). The numerical values for q ( y, t ) are found by negating thenumerically differentiated soliton centre location in x with respect to y .33igure 17: Full initial conditions for simulation of truncated (left) and bent-stemsolitons (right), satisfying the constraint (cid:82) u yy d x = 0. References [1] M. J. Ablowitz and D. E. Baldwin. Nonlinear shallow ocean-wave solitoninteractions on flat beaches.
Phys. Rev. E , 86(3):036305, September 2012.[2] M. J. Ablowitz, G. Biondini, and Q. Wang. Whitham modulation theory forthe Kadomtsev-Petviashvili equation.
Proc. R. Soc. A , 473(2204):20160695,2017.[3] M. J. Ablowitz and P. A. Clarkson.
Solitons, Nonlinear Evolution Equationsand Inverse Scattering . Cambridge University Press, 1991.[4] M. J. Ablowitz, A. Demirci, and Y.-P. Ma. Dispersive shock waves inthe Kadomtsev-Petviashvili and two dimensional Benjamin-Ono equations.
Physica D , 333:84–98, October 2016.[5] M. J. Ablowitz and H. Segur. On the evolution of packets of water waves.
J. Fluid Mech. , 92(04):691–715, 1979.[6] M. J. Ablowitz and J. Villarroel. On the Kadomtsev-Petviashvili Equationand Associated Constraints.
Stud. Appl. Math. , 85(3):195–213, 1991.[7] G. Biondini. Line Soliton Interactions of the Kadomtsev-Petviashvili Equa-tion.
Phys. Rev. Lett. , 99(6), August 2007.[8] G. Biondini, M. A. Hoefer, and A. Moro. Integrability, exact reductionsand special solutions of the KP-Whitham equations.
Nonlinearity , 33:4114–4132, 2020.[9] G. Biondini, K.-I. Maruno, M. Oikawa, and H. Tsuji. Soliton interactionsof the Kadomtsev-Petviashvili equation and generation of large-amplitudewater waves.
Stud. Appl. Math. , 122(4):377–394, 2009.3410] G. Biondini and T. Trogdon. Gibbs phenomenon for dispersive PDEs onthe line.
SIAM J. Appl. Math. , 77(3):813–837, 2017.[11] S. Chakravarty and Y. Kodama. Soliton solutions of the KP equation andapplication to shallow water waves.
Stud. Appl. Math. , 123(1):83–151, 2009.[12] S. Chakravarty, T. McDowell, and M. Osborne. Numerical studies of theKP line-solitons.
Comm. Nonlin. Sci. Num. Sim. , 44:37–51, March 2017.[13] R. Courant and K. O. Friedrichs.
Supersonic Flow and Shock Waves .Springer-Verlag, 1948.[14] M. Funakoshi. Reflection of obliquely incident solitary waves.
J. Phys. Soc.Japan , 49(6):2371–2379, 1980.[15] T. Grava, C. Klein, and G. Pitton. Numerical study of the Kadomtsev-Petviashvili equation and dispersive shock waves.
Proc. R. Soc. A ,474(2210):20170458, 2018.[16] R. Grimshaw. Evolution equations for long, nonlinear internal waves instratified shear flows.
Stud. Appl. Math. , 65(2):159–188, 1981.[17] C. R. Jackson. An atlas of internal solitary-like waves and their properties.Technical report, Global Ocean Associates, 2004.[18] B. B. Kadomtsev and V. I. Petviashvili. On the stability of solitary wavesin weakly dispersing media.
Sov. Phys. Dokl. , 15(6):753–756, 1970.[19] C.-Y. Kao and Y. Kodama. Numerical study of the KP equation for non-periodic waves.
Math. Comp. Sim. , 82(7):1185–1218, March 2012.[20] C. Klein, C. Sparber, and P. Markowich. Numerical study of oscillatoryregimes in the KadomtsevPetviashvili equation.
J. Nonlin. Sci. , 17(5):429–470, Jun 2007.[21] K. Ko and H. H. Kuehl. Cylindrical and spherical KortewegdeVries solitarywaves.
Phys. Fluids , 22(7):1343–1348, 1979.[22] Y. Kodama. KP solitons in shallow water.
J. Phys. A , 43(43):434004, 2010.[23] Y. Kodama, M. Oikawa, and H. Tsuji. Soliton solutions of the KP equationwith V-shape initial waves.
J. Phys. A , 42(31):312001, 2009.[24] Y. Kodama and H. Yeh. The KP theory and Mach reflection.
J. FluidMech. , 800:766–786, 2016.[25] Yuji Kodama. Young diagrams and N-soliton solutions of the KP equation.
J. Phys. A , 37(46):11169–11190, 2004.[26] P. Krehl and M. van der Geest. The discovery of the Mach reflection effectand its demonstration in an auditorium.
Shock Waves , 1:3–15, 1991.3527] S.-J. Lee and R. H. J. Grimshaw. Upstreamadvancing waves generated bythreedimensional moving disturbances.
Phys. Fluids A , 2(2):194–201, 1990.[28] W. Li, H. Yeh, and Y. Kodama. On the Mach reflection of a solitary wave:revisited.
J. Fluid Mech. , 672:326–357, 2011.[29] Y. Li and P. D. Sclavounos. Three-dimensional nonlinear solitary wavesin shallow water generated by an advancing disturbance.
J. Fluid Mech. ,470:383–410, 2002.[30] E. Mach and J. Wosyka. ¨Uber einige mechanische Wirkungen des elec-trischen Funkens.
Sitzungsber. Akad. Wiss. Wien , 72:44–52, 1875.[31] W. K. Melville. On the Mach reflexion of a solitary wave.
J. Fluid Mech. ,98(2):285–297, 1980.[32] J. W. Miles. Obliquely interacting solitary waves.
J. Fluid Mech. ,79(01):157–169, 1977.[33] J. W. Miles. Resonantly interacting solitary waves.
J. Fluid Mech. ,79(01):171–179, 1977.[34] A Nakamura and H.-H. Chen. Soliton solutions of the cylindrical KdVequation.
J. Phys. Soc. Japan , 50(2):711–718, 1981.[35] J. C. Neu.
Singular Perturbation in the Physical Sciences . American Math-ematical Society, Providence, RI, 2015.[36] J. Pan, D. A. Jay, and P. M. Orton. Analyses of internal solitary waves gen-erated at the Columbia River plume front using SAR imagery.
J. Geophys.Res. Oceans , 112(C7), 2007.[37] P. H. Perroud. The solitary wave reflection along a straight vertical wallat oblique incidence. Institute of Engineering Research. Wave ResearchLaboratory. Technical report ser. 99, University of California, Berkeley,1957.[38] J. Smoller.
Shock waves and reaction-diffusion equations . Springer, 1994.[39] M. Tanaka. Mach reflection of a large-amplitude solitary wave.
J. FluidMech. , 248:637–661, 1993.[40] H. Tsuji and M. Oikawa. Oblique interaction of solitons in an extendedKadomtsev-Petviashvili equation.
J. Phys. Soc. Japan , 76(8):084401, 2007.[41] C. Wang and R. Pawlowicz. Oblique wave-wave interactions of nonlin-ear near-surface internal waves in the Strait of Georgia.
J. Geophys. Res.Oceans , 117(C6), 2012.[42] C. Wang and R. Pawlowicz. Internal wave generation from tidal flow exitinga constricted opening.
J. Geophys. Res. Oceans , 122(1):110–125, 2017.3643] C. Yuan, R. Grimshaw, E. Johnson, and Z. Wang. Topographic effect onoblique internal wavewave interactions.