On Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction
OOn Coupled Boussinesq Equations and Ostrovsky equationsfree from zero-mass contradiction
K R Khusnutdinova and M R Tranter Department of Mathematical Sciences, Loughborough University, Loughborough LE11 3TU, UK Department of Mathematics and Physics, Nottingham Trent University, Nottingham NG11 8NS,UKE-mail:
Abstract.
Long weakly-nonlinear waves in a layered waveguide with an imperfect interface (softbonding between the layers) can be modelled using coupled Boussinesq equations. Previously, weconsidered the case when the materials of the layers have close mechanical properties, and the systemsupports radiating solitary waves. Here we are concerned with a more challenging case, when themechanical properties of the materials of the layers are significantly different, and the system supportswave packet solutions. We construct a weakly-nonlinear solution of the Cauchy problem for thissystem, considering the problem in the class of periodic functions on an interval of finite length. Thesolution is constructed for the deviation from the evolving mean value using a novel multiple-scalesprocedure involving fast characteristic variables and two slow time variables. By construction, theOstrovsky equations emerging within the scope of this procedure are solved for initial conditionswith zero mean while initial conditions for the original system may have non-zero mean values.Asymptotic validity of the solution is carefully examined numerically. We also discuss the applicationof the solution to the study of co-propagating waves generated by the solitary or cnoidal wave initialconditions, as well as the case of counter-propagating waves and the resulting wave interactions. Onelocal and two nonlocal conservation laws ae obtained and used to control the accuracy of the numericalsimulations.
Submitted to:
Nonlinearity a r X i v : . [ m a t h . A P ] F e b n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction The propagation of long weakly-nonlinear longitudinal bulk strain solitary waves in an elasticwaveguides can be modelled using Boussinesq-type equations (see [1, 2, 3, 4, 5, 6, 7, 8, 9]). CoupledBoussinesq equations were derived as a model describing similar waves in layered waveguides with asoft bonding layer [10]. The stability of solitary waves has inspired research into their applicabilityfor detection of delamination in some cases, in addition to existing methods [11, 12, 13, 14, 15].Low-frequency wave propagation in solids is relevant to a large number of modern applications (see,for example, [3, 4, 7, 16, 17, 18, 19] and references therein).Of interest in this paper is the system of coupled regularised Boussinesq (cRB) equations,presented here in non-dimensional and scaled form: u tt − u xx = ε (cid:20) (cid:0) u (cid:1) xx + u ttxx − δ ( u − w ) (cid:21) , (1.1) w tt − c w xx = ε (cid:104) α (cid:0) w (cid:1) xx + βw ttxx + γ ( u − w ) (cid:105) . (1.2)This system of equations has been derived to describe long nonlinear longitudinal bulk strain wavesin a bi-layer with a sufficiently soft bonding layer [10]. In this context, u and w denote longitudinalstrains in the layers, α , β , γ , δ are coefficients depending on the mechanical and geometrical propertiesof a waveguide, c is the ratio of the characteristic linear wave speeds in the layers and ε is a smallamplitude parameter. The detection of delamination regions was considered in [14], for strain solitarywaves incident on the waveguide.In all of these contexts, the natural initial conditions have non-zero mean, while the associateduni-directional equations emerging in the construction of weakly-nonlinear solutions are typicallysolved numerically using pseudo-spectral schemes and periodic boundary conditions. However, recentwork has shown that solving such equations on a periodic domain may require careful attention ifthe initial conditions have non-zero mean [20, 21]. The Boussinesq-Klein-Gordon equation can befound from the cRB equations by taking the limit γ → w = 0. The weakly-nonlinear solutionderived via traditional procedure leads at leading order to two uni-directional Ostrovsky equations[22] that necessarily require zero-mean initial conditions. The existence of such a formal constraintis known as the ‘zero-mass (or zero-mean) contradiction’ [21] . However, in [21] it was shown thatthe derivation procedure can be modified in order to develop a weakly-nonlinear solution of theBoussinesq-Klein-Gordon equation for a deviation from the evolving non-zero mean. The earlierresults in [20] were developed at the level of Fourier expansions in the spatial variable, while theprocedure suggested in [21] can be viewed as a nonlinear extension of the d’Alembert’s solution. TheOstrovsky equations emerging within the scope of this procedure are solved for zero-mean initialconditions by construction, and the zero-mean contradiction is avoided.Considering our cRB system of equations, in the case when u = w the system reduces to theregularised Boussinesq equation u tt − u xx = ε (cid:20) (cid:0) u (cid:1) xx + u ttxx (cid:21) , (1.3) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction c − O ( ε ), pure solitary waves are replaced with long-living radiating solitary waves[10, 23, 24]. When the materials are significantly different, characterised by c − O (1), the solitarywaves are replaced by wave packets [23] governed by the Ostrovsky equations [22]. Therefore, as wehave Ostrovsky equations and coupled Ostrovsky equations, which necessarily require zero-meaninitial conditions, we need to consider how a weakly-nonlinear solution can be constructed that takesaccount of this restriction.The first case, when c − O ( ε ), was constructed in [25] and is briefly overviewed at thebeginning of Section 2. We will re-examine this solution in comparison with the more complicatedcase when c − O (1), which is the main topic of the present paper. In this second case theleading-order linear wave operators of the two equations have different characteristic variables, whilethe previous case can be brought to the asymptotically equivalent form where they are the same.The rest of the paper is organised as follows. In Section 2 we construct a weakly-nonlinearsolution of the Cauchy problem for the cRB equations (1.1) - (1.2) in the case c − O (1) ona periodic domain, considering asymptotic multiple-scales expansions for the deviation from theoscillating mean values, and using two sets of fast characteristic variables and two slow time variables.The validity of the solutions is examined in Section 3 by comparing the constructed weakly-nonlinear solution with direct numerical simulations of the Cauchy problem, and we discuss bothcases c − O ( ε ) and c − O (1). In Section 4 we use both direct numerical simulations andthe constructed weakly-nonlinear solutions to study the behaviour of the waves generated by thesolitary and cnoidal wave initial conditions, co- and counter-propagating, for the cases when thecharacteristic speeds are either close or significantly different.In Section 5 we discuss the local and non-local conservation laws used to control the accuracyof numerical simulations and we conclude our studies in Section 6. We solve the equation system (1.1) - (1.2) on the periodic domain x ∈ [ − L, L ]. The initial-value(Cauchy) problem is considered, and the initial conditions are written as u ( x,
0) = F ( x ) , u t ( x,
0) = V ( x ) , (2.1) w ( x,
0) = F ( x ) , w t ( x,
0) = V ( x ) . (2.2)Firstly, we integrate (1.1) - (1.2) in x over the period 2 L to obtain evolution equations of the formd d t (cid:90) L − L u ( x, t ) d x + εδ (cid:90) L − L ( u ( x, t ) − w ( x, t )) = 0 , (2.3)d d t (cid:90) L − L w ( x, t ) d x − εγ (cid:90) L − L ( u ( x, t ) − w ( x, t )) = 0 . (2.4) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction u and w as u ( t ) := 12 L (cid:90) L − L u ( x, t ) d x, w ( t ) := 12 L (cid:90) L − L w ( x, t ) d x, (2.5)we can solve (2.3) - (2.4) to obtain u = d + δd cos ωt + d t + δd sin ωt, (2.6) w = d − γd cos ωt + d t − γd sin ωt. (2.7)Using the initial conditions (2.1), (2.2) we can determine the values of the coefficients as d = γF + δF δ + γ , d = F − F δ + γ , d = γV + δV ω ( δ + γ ) , d = V − V ω ( δ + γ ) , (2.8)where we introduce the notation ω = (cid:112) ε ( δ + γ ) and F i = (cid:90) L − L F i ( x ) d x, V i = (cid:90) L − L V i ( x ) d x, i = 1 , . (2.9)To simplify the problem we will consider initial conditions that satisfy the condition d = d = 0,that is 12 L (cid:90) L − L V i d x = 0 , i = 1 , . (2.10)This condition appears naturally in many physical applications, and we impose it here in order tosimplify our derivations. As was shown in [21] for the Boussinesq-Klein-Gordon equation (1.3), suchconditions can be relaxed.We subtract (2.6) from u and (2.7) from w to obtain an equation with zero mean value, so weintroduce ˜ u = u − u and ˜ w = w − w to obtain the problem for deviations˜ u tt − ˜ u xx = ε (cid:20) (cid:0) ˜ u (cid:1) xx + ( d + δd cos ωt ) ˜ u xx + ˜ u ttxx − δ (˜ u − ˜ w ) (cid:21) , (2.11)˜ w tt − c ˜ w xx = ε (cid:104) α (cid:0) ˜ w (cid:1) xx + α ( d − γd cos ωt ) ˜ w xx + β ˜ w ttxx + γ (˜ u − ˜ w ) (cid:105) . (2.12)Note that this problem has variable coefficients, and the traditional procedure used to derive uni-directional models of the Korteweg-de Vries/Ostrovsky type is no longer applicable. The initialconditions become ˜ u ( x,
0) = ˜ F ( x ) = F ( x ) − u, ˜ u t ( x,
0) = ˜ V ( x ) = V ( x ) , (2.13)˜ w ( x,
0) = ˜ F ( x ) = F ( x ) − w, ˜ w t ( x,
0) = ˜ V ( x ) = V ( x ) , (2.14)and, by construction, have zero mean value. We omit tildes in all functions to simplify notation insubsequent sections. n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction c − O ( ε )The first case when c − O ( ε ) was considered in [25], so we only summarise the results here. Inthis case the waves are resonant and an initial solitary wave solution in both layers will evolve intoa radiating solitary wave ; a solitary wave with a co-propagating one-sided oscillatory tail [10]. Wenote that we have ( c − /ε = O (1) and therefore we can rearrange (2.12) to give u tt − u xx = ε (cid:20) (cid:0) u (cid:1) xx + ( d + δd cos ωt ) u xx + βu ttxx − δ ( u − w ) (cid:21) , (2.15) w tt − w xx = ε (cid:20) α (cid:0) w (cid:1) xx + α (cid:18) d − γd cos ωt + c − ε (cid:19) w xx + βw ttxx + γ ( u − w ) (cid:21) . (2.16)We look for a weakly-nonlinear solution with two slow time scales [21], which takes the form u ( x, t ) = f − ( ξ − , τ, T ) + f +1 ( ξ + , τ, T ) + √ εP ( ξ − , ξ + , τ, T ) + εQ ( ξ − , ξ + , τ, T )+ ε R ( ξ − , ξ + , τ, T ) + ε S ( ξ − , ξ + , τ, T ) + O (cid:16) ε (cid:17) , (2.17) w ( x, t ) = f − ( ξ − , τ, T ) + f +2 ( ξ + , τ, T ) + √ εP ( ξ − , ξ + , τ, T ) + εQ ( ξ − , ξ + , τ, T )+ ε R ( ξ − , ξ + , τ, T ) + ε S ( ξ − , ξ + , τ, T ) + O (cid:16) ε (cid:17) , (2.18)with characteristic and slow time variables ξ ± = x ± t, τ = √ εt, T = εt. As we are considering the solution on the periodic domain, the functions u and w are 2 L -periodicfunctions in x . Therefore we require that f − , and f +1 , are periodic in ξ − and ξ + respectively, andthat all terms in the asymptotic expansions are products of the functions f ± , and their derivatives.This implies that all functions are periodic in ξ − , ξ + at fixed ξ + , ξ − .Furthermore, the functions f ± , have zero mean, satisfying12 L (cid:90) L − L f ± , d ξ ± = 0 , (2.19)and therefore all functions in the expansion have zero mean. These conditions are consistent withthose established in [21].We substitute (2.17) into (2.15) and (2.18) into (2.16), collecting terms at increasing powers of √ ε to determine expressions for all functions in the expansion. Matching up to powers of ε , we findthat u ( x, t ) = f − + f +1 + √ ε (cid:0) g − + g +1 (cid:1) + ε (cid:0) h − + h +1 + f c (cid:1) + O (cid:16) ε (cid:17) , (2.20) w ( x, t ) = f − + f +2 + √ ε (cid:0) g − + g +2 (cid:1) + ε (cid:0) h − + h +2 + f c (cid:1) + O (cid:16) ε (cid:17) , (2.21)where the functions f ± , are described by (cid:16) ∓ f ± T + f ± f ± ξ ± + d f ± ξ ± + f ± ξ ± ξ ± ξ ± (cid:17) ξ ± = δ (cid:0) f ± − f ± (cid:1) , (cid:18) ∓ f ± T + αf ± f ± ξ ± + (cid:18) αd + c − ε (cid:19) f ± ξ ± + βf ± ξ ± ξ ± ξ ± (cid:19) ξ ± = γ (cid:0) f ± − f ± (cid:1) , (2.22) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction g ± , at O ( √ ε ) can be expressed in terms of the leading order functions f ± , as g ± = ± d δ ω sin (˜ ωτ ) f ± ξ ± , g ± = ∓ αd γ ω sin (˜ ωτ ) f ± ξ ± . (2.23)At O ( ε ) we find h ± = ˜ ω ρ f ± − δ ( ρ + ρ )2 f ± − ˜ ω ρ f ± ξ ± ξ ± + φ ± ( ξ ± , T ) , (2.24) h ± t = − ˜ ω ρ f ± + γ ( ρ + ρ )2 f ± − ˜ ω ρ f ± ξ ± ξ ± + φ ± ( ξ ± , T ) , (2.25)and h c = − (cid:18) f ξ − (cid:90) ξ + − L f +1 ( σ ) d σ + 2 f − f +1 + f +1 ξ + (cid:90) ξ − − L f − ( σ ) d σ (cid:19) , (2.26) h c = − α (cid:18) f ξ − (cid:90) ξ + − L f +2 ( σ ) d σ + 2 f − f +2 + f +2 ξ + (cid:90) ξ − − L f − ( σ ) d σ (cid:19) , (2.27)where we introduced the coefficients ρ = − d δ ω cos(˜ ωτ ) , ρ = − αd δ ω cos(˜ ωτ ) . (2.28)The functions φ ± , are described by (cid:16) ∓ φ ± T + (cid:0) f ± φ ± (cid:1) ξ ± + d φ ± ξ ± + φ ± ξ ± ξ ± ξ ± (cid:17) ξ ± = δ (cid:0) φ ± − φ ± (cid:1) + f ± T T ∓ f ± ξ ± ξ ± ξ ± T + ˜ ω ˜ θ f ± ξ ± ξ ± − ˜ θ δ + αγ ) f ± ξ ± ξ ± − ˜ θ (cid:16) f ± ξ ± (cid:17) ξ ± ξ ± , (2.29)and (cid:18) ∓ φ ± T + α (cid:0) f ± φ ± (cid:1) ξ ± + αd φ ± ξ ± + c − ε φ ± ξ ± + βφ ± ξ ± ξ ± ξ ± (cid:19) ξ ± = γ (cid:0) φ ± − φ ± (cid:1) + f ± T T ∓ βf ± ξ ± ξ ± ξ ± T + ˜ ω ˜ θ f ± ξ ± ξ ± + ˜ θ α ( δ + αγ ) f ± ξ ± ξ ± − α ˜ θ (cid:16) f ± ξ ± (cid:17) ξ ± ξ ± , (2.30)where we have introduced the modified non-secular coefficients˜ θ = d δ ω , ˜ θ = αd γ ω . (2.31)These are complemented with initial conditions f ± | T =0 = 12 (cid:18) ˜ F ( x ± t ) ± (cid:90) x ± t − L ˜ V ( σ ) d σ (cid:19) , f ± | T =0 = 12 (cid:18) ˜ F ( x ± t ) ± (cid:90) x ± t − L ˜ V ( σ ) d σ (cid:19) . (2.32) φ ± = 12 (cid:18) J ∓ (cid:90) ξ ± − L K ( σ ) d σ (cid:19) , φ ± = 12 (cid:18) J ∓ (cid:90) ξ ± − L K ( σ ) d σ (cid:19) , (2.33) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction J = − ˜ ω ρ (cid:0) f − + f +1 (cid:1) + δ ( ρ + ρ )2 (cid:0) f − + f +2 (cid:1) + ˜ ω ρ (cid:16) f − ξ − ξ − + f +1 ξ + ξ + (cid:17) − h c ,K = f − T + f +1 T + ˜ ω ρ f ± (cid:16) f − ξ − − f +1 ξ + (cid:17) + δ ( ρ + ρ )2 (cid:16) f − ξ − − f +2 ξ + (cid:17) + ˜ ω ρ (cid:16) f − ξ − ξ − ξ − − f +1 ξ + ξ + ξ + (cid:17) ,J = ˜ ω ρ (cid:0) f − + f +2 (cid:1) − γ ( ρ + ρ )2 (cid:0) f − + f +1 (cid:1) + ˜ ω ρ (cid:16) f − ξ − ξ − + f +2 ξ + ξ + (cid:17) − h c ,K = f − T + f +2 T − ˜ ω ρ (cid:16) f − ξ − − f +2 ξ + (cid:17) − γ ( ρ + ρ )2 f ± (cid:16) f − ξ − − f +1 ξ + (cid:17) + ˜ ω ρ (cid:16) f − ξ − ξ − ξ − − f +2 ξ + ξ + ξ + (cid:17) . (2.34) c − O (1)The focus of this paper is the case when the phase speeds are distinctly different i.e. c − O (1).In this case the characteristic variables cannot be the same in each layer, and instead we have twodistinct pairs of characteristic variables. Therefore we look for a weakly-nonlinear solution of theform u ( x, t ) = f − ( ξ − , τ, T ) + f +1 ( ξ + , τ, T ) + √ εP ( ξ − , ξ + , τ, T ) + εQ ( ξ − , ξ + , τ, T )+ ε R ( ξ − , ξ + , τ, T ) + ε S ( ξ − , ξ + , τ, T ) + O (cid:16) ε (cid:17) , (2.35) w ( x, t ) = f − ( ν − , τ, T ) + f +2 ( ν + , τ, T ) + √ εP ( ν − , ν + , τ, T ) + εQ ( ν − , ν + , τ, T )+ ε R ( ν − , ν + , τ, T ) + ε S ( ν − , ν + , τ, T ) + O (cid:16) ε (cid:17) , (2.36)where we have the following characteristic and slow time variables ξ ± = x ± t, ν ± = x ± ct, τ = √ εt, T = εt. We follow the same steps as for the previous case, by substituting (2.35) and (2.36) into (2.11) and(2.12) and comparing at increasing powers of √ ε , using the same assumptions of periodicity and zeromean for all functions in the expansion. The equations are satisfied at leading order so we move ontoterms at O ( √ ε ). At this order we have − P ξ − ξ + − f − ξ − τ + 2 f +1 ξ + τ = 0 . (2.37)We average (2.37) with respect to the fast spatial variable x at constant ξ − and ξ + (see [21]).Averaging P ξ − ξ + at constant ξ − gives12 L (cid:90) L − L P ξ − ξ + d x = 14 L (cid:90) L − ξ − − L − ξ − P ξ − ξ + d ξ + = 14 L (cid:2) P ξ − (cid:3) L − ξ − − L − ξ − = 0 , (2.38)with a similar result for averaging at constant ξ + . Applying the averaging to (2.37) we find f − ξ − τ = 0 = ⇒ f − = ˜ f − ( ξ − , T ) and f +1 ξ + τ = 0 = ⇒ f +1 = ˜ f +1 ( ξ + , T ) . (2.39) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction w we obtain − c P ν − ν + − cf − ν − τ + 2 cf +2 ν + τ = 0 , (2.40)which, after averaging, leads to f − = ˜ f − ( ν − , T ) and f +2 = ˜ f +2 ( ν + , T ) . (2.41)Substituting (2.39) into (2.37) we obtain P ξ − ξ + = 0 = ⇒ P = g − ( ξ − , τ, T ) + g +1 ( ξ + , τ, T ) , (2.42)and similarly substituting (2.41) into (2.40) we obtain P ν − ν + = 0 = ⇒ P = g − ( ν − , τ, T ) + g +2 ( ν + , τ, T ) . (2.43)As before, we omit the tildes on f ± , in subsequent steps. The initial condition for f ± , is found inthe same way as before, by substituting (2.35) into (2.1) and (2.36) into (2.2) and comparing termsat O (1) to obtain f − + f +1 (cid:12)(cid:12) T =0 = F ( x ) , − f − ξ − + f +1 ξ + (cid:12)(cid:12)(cid:12) T =0 = V ( x ) , = ⇒ f ± | T =0 = 12 (cid:18) F ( x ± t ) ± (cid:90) x ± t − L V ( σ ) d σ (cid:19) , (2.44)and (cid:40) f − + f +2 (cid:12)(cid:12) T =0 = F ( x ) , − cf − ν − + cf +2 ν + (cid:12)(cid:12) T =0 = V ( x ) , = ⇒ f ± | T =0 = 12 c (cid:18) cF ( x ± ct ) ± (cid:90) x ± ct − L V ( σ ) d σ (cid:19) . (2.45)Comparing terms at O ( ε ), using the results from the previous order, we have − Q ξ − ξ + = 2 g − ξ − τ + (cid:16) f − T + f − f − ξ − + d f − ξ − + f − ξ − ξ − ξ − (cid:17) ξ − − δ (cid:0) f − − f − (cid:1) − g +1 ξ + τ + (cid:16) − f +1 T + f +1 f +1 ξ + + d f +1 ξ + + f +1 ξ + ξ + ξ + (cid:17) ξ + − δ (cid:0) f +1 − f +2 (cid:1) + d δ cos (˜ ωτ ) (cid:16) f − ξ − ξ − + f +1 ξ + ξ + (cid:17) + f − ξ − ξ − f +1 + 2 f − ξ − f +1 ξ + + f − f +1 ξ + ξ + , (2.46)for the equation governing u , and for the equation governing w we obtain − c Q ν − ν + = 2 cg − ν − τ + (cid:0) cf − T + αf − f − ν − + αd f − ν − + βc f − ν − ν − ν − (cid:1) ν − + γ (cid:0) f − − f − (cid:1) − cg +2 ν + τ + (cid:0) − cf +2 T + αf +2 f +2 ν + + αd f +2 ν + + βc f +2 ν + ν + ν + (cid:1) ν + + γ (cid:0) f +1 − f +2 (cid:1) − αd γ cos (˜ ωτ ) (cid:0) f − ν − ν − + f +2 ν + ν + (cid:1) + α (cid:0) f − ν − ν − f +2 + 2 f − ν − f +2 ν + + f − f +2 ν + ν + (cid:1) . (2.47)Averaging (2.46) at constant ξ − or constant ξ + gives ± g ξ ± τ = d δ cos (˜ ωτ ) f ± ξ ± ξ ± + A ( ξ ± , T ) , (2.48) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction A = (cid:16) ∓ f ± T + f ± f ± ξ ± + d f ± ξ ± + f ± ξ ± ξ ± ξ ± (cid:17) ξ ± − δf ± . (2.49)To avoid secular terms we require that A = 0. Therefore we obtain an Ostrovsky equation for f ± of the form (cid:16) ∓ f ± T + f ± f ± ξ ± + d f ± ξ ± + f ± ξ ± ξ ± ξ ± (cid:17) ξ ± = δf ± . (2.50)Similarly, averaging (2.47) at constant ν − or constant ν + leads to the equation ± cg ν ± τ = − d γ cos (˜ ωτ ) f ± ν ± ν ± + A ( ν ± , T ) , (2.51)where A = (cid:0) ∓ cf ± T + αf ± f ± ν ± + αd f ± ν ± + βc f ± ν ± ν ± ν ± (cid:1) ν ± − γf ± . (2.52)We again require A = 0 to avoid secular terms, and therefore we obtain an Ostrovsky equation for f ± of the form (cid:0) ∓ cf ± T + αf ± f ± ν ± + αd f ± ν ± + βc f ± ν ± ν ± ν ± (cid:1) ν ± = γf ± . (2.53)As before we integrate (2.48) and (2.51), taking account of (2.50) and (2.53), to find an equation for g ± , of the form g ± = ± d δ ω sin (˜ ωτ ) f ± ξ ± + G ± ( ξ ± , T ) = ± θ f ± ξ ± + G ± ( ξ ± , T ) , (2.54)and g ± = ∓ αd γ c ˜ ω sin (˜ ωτ ) f ± ν ± + G ± ( ν ± , T ) = ∓ θ f ± ν ± + G ± ( ν ± , T ) , (2.55)where G ± and G ± are functions to be found, and we introduce θ = d δ ω sin (˜ ωτ ) , θ = αd γ c ˜ ω sin (˜ ωτ ) . (2.56)Substituting (2.50) and (2.54) into (2.46) and integrating we obtain Q = h − ( ξ − , τ, T ) + h +1 ( ξ + , τ, T ) + h c ( ξ − , ξ + , T ) + ˆ f − ( ν − , T ) + ˆ f +2 ( ν + , T ) , (2.57)where h c = − (cid:18) f ξ − (cid:90) ξ + − L f +1 ( σ ) d σ + 2 f − f +1 + f +1 ξ + (cid:90) ξ − − L f − ( σ ) d σ (cid:19) , (2.58)and ˆ f ± = δc − (cid:90) ν ± − L (cid:90) v − L f ± ( u, T ) d u d v. (2.59)Similarly substituting (2.53) and (2.55) into (2.47) and integrating we obtain Q = h − ( ν − , τ, T ) + h +2 ( ν + , τ, T ) + h c ( ν − , ν + , T ) + ˆ f − ( ξ − , T ) + ˆ f +1 ( ξ + , T ) , (2.60)where h c = − α c (cid:18) f ν − (cid:90) ν + − L f +2 ( σ ) d σ + 2 f − f +2 + f +2 ν + (cid:90) ν − − L f − ( σ ) d σ (cid:19) , (2.61) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction f ± = − γc − (cid:90) ξ ± − L (cid:90) v − L f ± ( u, T ) d u d v. (2.62)As before we find the initial condition for G ± , by substituting (2.35) and (2.36) into (2.1) and (2.2),comparing terms at O ( √ ε ). Therefore, taking account of the results found for (2.54) and (2.55), andnoting that θ | T =0 = 0, we obtain θ f − ξ − + θ f +1 ξ + + G − + G +1 (cid:12)(cid:12)(cid:12) T =0 = 0 , − θ f − ξ − ξ − + θ f +1 ξ + ξ + − G − ξ − + G +1 ξ + (cid:12)(cid:12)(cid:12) T =0 = 0 , = ⇒ G ± (cid:12)(cid:12) T =0 = 0 . (2.63)Similarly we have (cid:40) θ f − ν − + θ f +2 ν + + G − + G +2 (cid:12)(cid:12) T =0 = 0 , − cθ f − ν − ν − + cθ f +2 ν + ν + − cG − ν − + cG +2 ν + (cid:12)(cid:12) T =0 = 0 , = ⇒ G ± (cid:12)(cid:12) T =0 = 0 . (2.64)We now aim to find equations for h ± , . Comparing terms at O (cid:0) ε / (cid:1) in the expansion of (2.11) and(2.12), taking account of results at previous orders, we have − R ξ − ξ + = 2 h − ξ − τ − h ξ + τ + (cid:16) g − T + (cid:0) f − g − (cid:1) ξ − + d g − ξ − + g − ξ − ξ − ξ − (cid:17) ξ − − δ (cid:0) g − − g − (cid:1) − g − ττ − g +1 ττ + (cid:16) − g +1 T + (cid:0) f +1 g +1 (cid:1) ξ + + d g +1 ξ + + g +1 ξ + ξ + ξ + (cid:17) ξ + − δ (cid:0) g +1 − g +2 (cid:1) + d δ cos (˜ ωτ ) (cid:16) g − ξ − ξ − + g +1 ξ + ξ + (cid:17) + g − ξ − ξ − f +1 + 2 g − ξ − f +1 ξ + + g − f +1 ξ + ξ + + g +1 ξ + ξ + f − + 2 g +1 ξ + f − ξ − + g +1 f − ξ − ξ − , (2.65)and − c R ν − ν + = 2 ch − ν − τ + (cid:16) cg − T + α (cid:0) f − g − (cid:1) ν − + αd g − ν − + βc g − ν − ν − ν − (cid:17) ν − + γ (cid:0) g − − g − (cid:1) − ch +2 ν + τ + (cid:16) − cg +2 T + α (cid:0) f +2 g +2 (cid:1) ν + + αd g +2 ν + + βc g +2 ν + ν + ν + (cid:17) ν + + γ (cid:0) g +1 − g +2 (cid:1) − αd γ cos (˜ ωτ ) (cid:0) g − ν − ν − + g +2 ν + ν + (cid:1) − g − ττ − g +2 ττ + α (cid:2) g − ν − ν − f +2 + 2 g − ν − f +2 ν + + g − f +2 ν + ν + + g +2 ν + ν + f − + 2 g +2 ν + f − ν − + g +2 f − ν − ν − (cid:3) . (2.66)Substituting (2.54) into (2.65) and averaging at constant ξ − or constant ξ + leads to ± h ± ξ ± τ = ± θ (cid:16) ∓ f ± T + f ± f ± ξ ± + d f ± ξ ± + f ξ ± ξ ± ξ ± (cid:17) ξ ± ξ ± ∓ θ δf ± + (cid:16) ∓ G ± T + (cid:0) f ± G ± (cid:1) ξ ± + d G ± ξ ± + G ξ ± ξ ± ξ ± (cid:17) ξ ± − δG ± ξ ± ± θ ˜ ω f ± ξ ± ± θ d δ cos (˜ ωτ ) f ± ξ ± ξ ± ξ ± . (2.67)If we differentiate (2.50) with respect to the appropriate characteristic variable, we can eliminate thefirst line from (2.67) to obtain an expression for h ± ξ ± τ of the form2 h ± ξ ± τ = θ ˜ ω f ± ξ ± + θ d δ cos (˜ ωτ ) f ± ξ ± ξ ± ξ ± + ˜ G ± ( ξ ± , T ) , (2.68) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction G ± ( ξ ± , T ) = (cid:16) ∓ G ± T + (cid:0) f ± G ± (cid:1) ξ ± + d G ± ξ ± + G ± ξ ± ξ ± ξ ± (cid:17) ξ ± − δG ± . (2.69)To avoid secular terms we require that ˜ G ± = 0 and therefore we have an equation for G ± of the form (cid:16) ∓ G ± T + (cid:0) f ± G ± (cid:1) ξ ± + d G ± ξ ± + G ± ξ ± ξ ± ξ ± (cid:17) ξ ± = δG ± . (2.70)Taking account of the initial condition in (2.63) and the form of (2.70), we can clearly see that G ± ≡ G ± . Integrating(2.68) we obtain h ± = − ˜ ω ρ f ± − ˜ ω ρ f ± ξ ± ξ ± + φ ± ( ξ ± , T ) , (2.71)where the functions φ ± are to be found at the next order and ρ = ∂ − θ = d δ ω cos (˜ ωτ ) . (2.72)Averaging (2.66) at constant ν − or ν + , and using (2.55) and (2.53) as was done above, we obtain anexpression for h ± of the form ± ch ± ξ ± τ = ∓ θ ˜ ω f ± ν ± ± θ d γ cos (˜ ωτ ) f ± ν ± ν ± ν ± + ˜ G ± ( ν ± , T ) , (2.73)where ˜ G ± ( ξ ± , T ) = (cid:16) ∓ cG ± T + α (cid:0) f ± G ± (cid:1) ν ± + αd G ± ν ± + βc G ± ν ± ν ± ν ± (cid:17) ν ± − γG ± . (2.74)As before we require that ˜ G ± = 0 and G ± is governed by the equation (cid:16) ∓ cG ± T + α (cid:0) f ± G ± (cid:1) ν ± + αd G ± ν ± + βc G ± ν ± ν ± ν ± (cid:17) ν ± = γG ± . (2.75)By the same argument as for G ± we have G ± ≡ h ± = ˜ ω ρ c f ± − ˜ ω ρ f ± ν ± ν ± + φ ± ( ν ± , T ) , (2.76)where again we need to find the function φ ± and ρ = ∂ − θ = αd γ ω cos (˜ ωτ ) . (2.77)Substituting (2.71) into (2.65) and integrating with respect to the appropriate characteristic variableswe find R = ψ − ( ξ − , τ, T ) + ψ +1 ( ξ + , τ, T ) + ψ c ( ξ − , ξ + , T ) + ˆ g − + ˆ g +2 , (2.78)where ψ c = − θ (cid:20) f +1 ξ + ξ + (cid:90) ξ − − L f − ( σ ) d σ − f +1 f − ξ − + f − f +1 ξ + − f − ξ − ξ − (cid:90) ξ + − L f +1 ( σ ) d σ (cid:21) , (2.79) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction g ± can be found by replacing f with g in (2.59). In a similar way we find an expression for R of the form R = ψ − ( ν − , τ, T ) + ψ +2 ( ν + , τ, T ) + ψ c ( ξ − , ξ + , T ) + ˆ g − + ˆ g +1 , (2.80)where ψ c = − αθ c (cid:20) f − ξ − ξ − (cid:90) ξ + − L f +2 ( σ ) d σ − f − ξ − f +2 + f − f +2 ξ + − f +2 ξ + ξ + (cid:90) ξ − − L f − ( σ ) d σ (cid:21) , (2.81)and again the expression for ˆ g ± can be found by replacing f with g in (2.62). The initial conditionfor the function φ ± is found by again substituting (2.35) into (2.1) and comparing terms at O ( ε ),taking account of (2.71). Therefore we obtain for φ ± h − + h +1 + h c + ˆ f − + ˆ f +2 (cid:12)(cid:12)(cid:12) T =0 = 0 ,f − T + f +1 T + g − τ + g +1 τ − h − ξ − + h +1 ξ + − h cξ − + h cξ + − c ˆ f − ξ − + c ˆ f +2 ξ + (cid:12)(cid:12)(cid:12) T =0 = 0 , = ⇒ φ ± = 12 (cid:18) J ± ∓ (cid:90) ξ ± − L K ( σ ) d σ (cid:19) , (2.82)where J = ˜ ω ρ (cid:0) f − + f +1 (cid:1) + ˜ ω ρ (cid:16) f − ξ − ξ − + f +1 ξ + ξ + (cid:17) − h c − δc − (cid:90) ν − − L (cid:90) v − L f − ( u, T ) d u d v − δc − (cid:90) ν + − L (cid:90) v − L f +2 ( u, T ) d u d v,K = f − T + f +1 T − ˜ ω ρ (cid:16) f − ξ − − f +1 ξ + (cid:17) + ˜ ω ρ (cid:16) f − ξ − ξ − ξ − − f +1 ξ + ξ + ξ + (cid:17) − h cξ − + h cξ + − cδc − (cid:90) ν − − L f − ( u, T ) d u + cδc − (cid:90) ν + − L f +2 ( u, T ) d u. (2.83)Similarly for φ ± we find the initial condition by substituting (2.36) into (2.2) and comparing termsat O ( ε ), using (2.76). This gives h − + h +2 + h c + ˆ f − + ˆ f +1 (cid:12)(cid:12)(cid:12) T =0 = 0 ,f − T + f +2 T + g − τ + g +2 τ − ch − ν − + ch +2 ν + − ch cν − + ch cν + − ˆ f − ν − + ˆ f +1 ν + (cid:12)(cid:12)(cid:12) T =0 = 0 , = ⇒ φ ± = 12 c (cid:18) J ± ∓ (cid:90) ξ ± − L K ( σ ) d σ (cid:19) , (2.84)where J = − ˜ ω ρ c (cid:0) f − + f +2 (cid:1) + ˜ ω ρ (cid:0) f − ν − ν − + f +2 ν + ν + (cid:1) − h c + γc − (cid:90) ξ − − L (cid:90) v − L f − ( u, T ) d u d v + γc − (cid:90) ξ + − L (cid:90) v − L f +1 ( u, T ) d u d v,K = f − T + f +2 T + ˜ ω ρ (cid:0) f − ν − − f +2 ν + (cid:1) + c ˜ ω ρ (cid:0) f − ν − ν − ν − − f +2 ν + ν + ν + (cid:1) − ch cν − + ch cν + + γc − (cid:90) ξ − − L f − ( u, T ) d u − γc − (cid:90) ξ + − L f +1 ( u, T ) d u. (2.85) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction O ( ε ), however we still need to find an equation governing φ ± , ,therefore we compare terms at O ( ε ). The weakly-nonlinear expansions (2.35) and (2.36) aresubstituted into (2.15) and (2.16) respectively, taking account of all previous results. All couplingterms between left- and right-propagating functions are gathered in one function, and terms of thetype of (2.59) and (2.62) are gathered in another function for convenience, as we do not require themto determine φ ± , . Gathering terms at O ( ε ) we have − S ξ − ξ + = − f − T T − f +1 T T − g − τT − g +1 τT − h − ττ − h +1 ττ + 2 h − ξ − T − h +1 ξ + T + 2 ψ − ξ − τ − ψ +1 ξ + τ + (cid:0) f − h − (cid:1) ξ − ξ − + (cid:0) f +1 h +1 (cid:1) ξ + ξ + + 12 (cid:16) g − (cid:17) ξ − ξ − + 12 (cid:16) g + (cid:17) ξ + ξ + + h − ξ − ξ − ξ − ξ − + h +1 ξ + ξ + ξ + ξ + + ( d + d δ cos (˜ ωτ )) (cid:16) h − ξ − ξ − + h +1 ξ + ξ + (cid:17) − g − ξ − ξ − ξ − τ + 2 g +1 ξ + ξ + ξ + τ − f − ξ − ξ − ξ − T + 2 f +1 ξ + ξ + ξ + T − δ (cid:0) h − − h − + h +1 − h +2 (cid:1) − µ c − , (2.86)where µ c is the coupling terms at this order and Υ is the terms involving ˆ f or equivalent. For thesecond equation we have − c S ν − ν + = − f − T T − f +2 T T − g − τT − g +2 τT − h − ττ − h +2 ττ + 2 ch − ξ − T − ch +2 ξ + T + 2 cψ − ξ − τ − cψ +2 ξ + τ + α (cid:0) f − h − (cid:1) ξ − ξ − + α (cid:0) f +2 h +2 (cid:1) ξ + ξ + + α (cid:16) g − (cid:17) ξ − ξ − + α (cid:16) g + (cid:17) ξ + ξ + + ( αd − αd γ cos (˜ ωτ )) (cid:16) h − ξ − ξ − + h +2 ξ + ξ + (cid:17) + βc h − ξ − ξ − ξ − ξ − + βc h +2 ξ + ξ + ξ + ξ + − βcg − ξ − ξ − ξ − τ + 2 βcg +2 ξ + ξ + ξ + τ − βcf − ξ − ξ − ξ − T + 2 βcf +2 ξ + ξ + ξ + T + γ (cid:0) h − − h − + h +1 − h +2 (cid:1) − µ c − , (2.87)where again µ c is the coupling terms at this order of the expansion and Υ is the terms involvingˆ f or equivalent. As was done for the first case, we average (2.86) at constant ξ − or constant ξ + , oraverage (2.87) at constant ν − or constant ν + , integrate with respect to τ and rearrange to obtain ± ψ ξ ± = H ± ( ξ ± , τ, T ) + ˆ H ± ( ξ ± , T ) τ and ± cψ ν ± = H ± ( ν ± , τ, T ) + ˆ H ± ( ν ± , T ) τ, (2.88)where the functions H ± , , ˆ H ± , can be found from (2.86) and (2.87). To avoid secular terms we requirethat ˆ H ± , = 0 and this allows us to find equations for φ ± , . Therefore we look for terms in (2.86)that depend only on ξ ± and T . Similarly, we look for terms in (2.87) that only depend on ν ± and T .Following this approach we obtain the equations (cid:16) ∓ φ ± T + (cid:0) f ± φ ± (cid:1) ξ ± + d φ ± ξ ± + φ ± ξ ± ξ ± ξ ± (cid:17) ξ ± = δφ ± + f ± T T ∓ f ± ξ ± ξ ± ξ ± T (2.89)+ ˜ ω ˜ θ f ± ξ ± ξ ± − ˜ θ (cid:16) f ± ξ ± (cid:17) ξ ± ξ ± , (2.90)and (cid:16) ∓ cφ ± T + α (cid:0) f ± φ ± (cid:1) ξ ± + αd φ ± ξ ± + βc φ ± ξ ± ξ ± ξ ± (cid:17) ξ ± = γφ ± + f ± T T ∓ cβf ± ξ ± ξ ± ξ ± T + ˜ ω ˜ θ f ± ξ ± ξ ± − α ˜ θ (cid:16) f ± ξ ± (cid:17) ξ ± ξ ± , (2.91) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction θ = θ sin (˜ ωτ ) = d δ ω , ˜ θ = θ sin (˜ ωτ ) = αd γ c ˜ ω . (2.92)We have now defined all functions up to and including O ( ε ) and so stop our derivation, however intheory this could be continued to any order. The weakly-nonlinear solutions presented in Sections 2.1 and 2.2 are different, despite the governingequations being the same. Namely, we see the following key differences: • In the case when c − O ( ε ) we have just one set of characteristic variables, while in the casewhen c − O (1) we have to develop our asymptotic expansions using two distinct sets ofcharacteristic variables. • In the case when c − O ( ε ) we derive coupled Ostrovsky equations to describe the evolution ofthe functions f ± , , while in the case when c − O (1) we have uncoupled Ostrovsky equations. • Similarly, for the functions φ ± , , we derive linearised coupled Ostrovsky equations for c − O ( ε )and linearised uncoupled Ostrovsky equations for the case of c − O (1). • In the terms at O ( ε ), there is a coupling term in both cases. For the case when c − O ( ε ) wehave a linear term in f ± for u and vice-versa for w . However, for the case when c − O (1),we have a double integral term instead. • There is an additional term in the right-hand side of the equation for φ ± , in the case when c − O ( ε ), as the characteristic variables are the same so it survives the averaging.Importantly, in both cases the constructed solutions are free from the zero-mass contradictionsince the Ostrovsky equations have been derived for the zero-mean functions by construction. In Section 2 we constructed the weakly-nonlinear solution for the case of close and distinctcharacteristic speeds. We now confirm the validity of the constructed expansions by numericallysolving the system (1.1) - (1.2) and comparing this direct numerical solution to the constructedsolution (2.17) and (2.18) with an increasing number of terms included. The first case, when c − O ( ε ), was analysed in [25]. Here we determine the validity of the expansion in the secondcase, when c − O (1). This was constructed in Section 2.2, so we need to solve (2.39), (2.41) forthe leading-order solution and (2.90), (2.91) for the solution up to and including terms at O ( ε ).In order to numerically solve these equations, in this section and the subsequent section, weimplement three pseudospectral numerical schemes. For the coupled Boussinesq equations we usethe technique outlined in [25, 26], while for the case of a single Ostrovsky equation we use the methodin [21]. In [25] a pseudospectral method was outlined for coupled Ostrovsky equations, however wepresent a modified scheme here, based upon [27], to allow for larger time steps to be taken. Thisis presented in Appendix A. Unless stated otherwise, we assume that ∆ x = 0 .
1, ∆ t = 0 .
01 and n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction T = ε ∆ t . In some calculations these parameters may be changed to obtain a higher accuracyresult and this will be stated in the figure captions. The domain is taken as [ − L, L ] in all cases, with L = 300 for most calculations with solitary waves and 3 × L K for cnoidal waves, where L K is theperiod of the cnoidal wave.In all subsequent calculations, we shall refer to the solution of the system (1.1) - (1.2) as the“exact solution” and compare this to the weakly-nonlinear solution with an increasing number ofterms included. We calculate the solution in the domain x ∈ [ − ,
40] and for t ∈ [0 , T ] where T = 1 /ε . The initial conditions are taken to be F ( x ) = A sech (cid:18) x Λ (cid:19) + p, F ( x ) = A sech (cid:18) x Λ (cid:19) ,V ( x ) = 2 A Λ1 sech (cid:18) x Λ (cid:19) tanh (cid:18) x Λ (cid:19) , V ( x ) = 2 c A Λ sech (cid:18) x Λ (cid:19) tanh (cid:18) x Λ (cid:19) , (3.1)where p is a constant and we have A = 6 k , Λ = √ /k , k = 1 / √ A = 6 ck /α , Λ = √ cβ/k , k = (cid:112) α/ c . We have only added a pedestal to the initial condition for u in view of the translationsymmetry of the system. This also gives distinct non-zero values for d and d .We choose α = β = c = 2, δ = γ = 0 . p = 7 and present the comparison between the exactand weakly-nonlinear solutions at various orders of ε in Figure 1. We see from the enhanced image -40 -20 0 20 40456 -3.5 -3 -2.53.93.95 (a) γ = 0 . p = 7 for u . -40 -20 0 20 40123 -9 -8 -71.121.141.16 (b) γ = 0 . p = 7 for w . Figure 1.
A comparison of the direct numerical simulation (solid, blue) and the weakly-nonlinear solution including leading order (dashed, red), O ( √ ε ) (dash-dot, black) and O ( ε )(dash-dot, green) corrections, at t = 1 /ε , for (a) u and (b) w . Parameters are L = 40, N = 800, k = 1 / √ α = β = c = 2, γ = 0 . ε = 0 . t = 0 .
01 and ∆ T = ε ∆ t .The solution agrees well to leading order and this agreement is improved with the additionof higher-order corrections. that the leading order solution (red, dashed line) is improved with the addition of the O ( √ ε ) terms(black, dash-dotted line), and further improved by the inclusion of O ( ε ) terms (green, dash-dottedline). We only present one example here for brevity, but the same results have been observed formultiple values of δ , γ and pedestal height p . n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction ε to show that the errors behave as expected. We denotethe direct numerical solution to the system (1.1) - (1.2) as u num , the weakly-nonlinear solution (2.35),(2.36) with only leading order terms included as u , with terms up to and including O ( √ ε ) as u and with terms up to and including O ( ε ) as u . The maximum absolute error over x is defined as e i = max − L ≤ x ≤ L | u num ( x, t ) − u i ( x, t ) | , i = 1 , , . (3.2)This error is calculated at every time step and, to smooth the oscillations in the errors we average the e i values in the final third of the calculation, denoting this value as ˆ e i . We then use a least-squarespower fit to determine how the maximum absolute error varies with the small parameter ε . Thereforewe write the errors in the form exp [ˆ e i ] = C i ε α i , (3.3)and take the logarithm of both sides to form the error plot (the exponential factor is included so thatwe have ˆ e i as the plotting variable).The corresponding errors for the cases considered in Figure 1 are plotted in Figure 2. The slope -8 -7 -6-10-5 (a) γ = 0 . p = 7 for u . -8 -7 -6-10-5 (b) γ = 0 . p = 7 for w . Figure 2.
A comparison of error curves for varying values of ε , at t = 1 /ε , for the weakly-nonlinear solution including leading order (upper, blue), O ( √ ε ) (middle, red) and O ( ε )(lower, black) corrections, for (a) u and (b) w . Parameters are L = 40, N = 800, k = 1 / √ α = β = c = 2, γ = 0 .
5, ∆ t = 0 .
01 and ∆ T = ε ∆ t . The inclusion of more terms in theexpansion increases the accuracy. of the curves is approximately 0.77, 0.99 and 1.96 for u and 0.77, 1.00 and 1.95 for w , which areslightly higher than the theoretical values for the inclusion of leading order and the case with terms upto and including O ( ε ) terms. However this can be explained by the solutions in this case being wavepackets, so a phase shift at the level of O ( √ ε ) will not have as much as an effect on the errors as itwould in the first case, when the solution was close to a solitary wave [25]. Furthermore, the increasein the slope is consistent with the results seen for previous studies for the Boussinesq-Klein-Gordonequation [21]. n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction We now consider various cases for wave propagation, comparing the numerical solution to the cRBequations to the constructed weakly-nonlinear solution with the inclusion of terms up to O ( √ ε ) totake account of the mass. Two types of initial condition are considered: solitary wave and cnoidalwave. In each case, the first example will be for waves propagating in the same direction (co-propagating waves) and the second example will show waves propagating in opposite directions andthen interacting with each other (counter-propagating waves). These results are presented for thecase when the characteristic speeds are close ( c − O ( ε )) and when the characteristic speeds aresignificantly different ( c − O (1)). Firstly we take solitary wave initial conditions, as was done in Section 3. These are presented in(3.1) but are reproduced here, where the pedestal term p is removed so that both waves are on azero pedestal. Namely, for the co-propagating case, we have F ( x ) = A sech (cid:18) x + x Λ (cid:19) , F ( x ) = A sech (cid:18) x + x Λ (cid:19) ,V ( x ) = 2 A Λ1 sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) , V ( x ) = 2 c A Λ sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) , (4.1)where A = 6 k , Λ = √ /k , A = 6 ck /α , Λ = √ cβ/k . We will define k and k for eachsimulation as they vary between cases. For the counter-propagating waves, we introduce a secondwave, with the same parameters as the first wave, at a phase shift with a sign change in V to resultin wave propagation to the left for the second wave. Explicitly we have F ( x ) = A sech (cid:18) x + x Λ (cid:19) + A sech (cid:18) x + x Λ (cid:19) , F ( x ) = A sech (cid:18) x + x Λ (cid:19) + A sech (cid:18) x + x Λ (cid:19) ,V ( x ) = 2 A Λ1 sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) − A Λ1 sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) ,V ( x ) = 2 c A Λ sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) − c A Λ sech (cid:18) x + x Λ (cid:19) tanh (cid:18) x + x Λ (cid:19) . (4.2) We consider the case when the characteristic speeds in theequations are close, so we have c − O ( ε ). The case of co-propagating waves was exploredin [25], so we exclude these from our considerations here. Instead we consider the case of a counter-propagating wave. We take the initial condition (4.2), where the waves will be well separated viachoice of x and x . The results are shown in Figure 3, where the evolution of the wave is shownfor u and w in a single plot. Firstly we see that the solitons form a co-propagating radiating tail,as expected from the co-propagating case [25]. As they interact the solitons appear to emergewith only a small change in their amplitudes, and there is good agreement between the direct and n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -300 -200 -100 0 100 200 30000.511.5 (a) Solution for u at t = 150, t = 250 and t = 500. -300 -200 -100 0 100 200 30000.511.5 (b) Solution for w at t = 150, t = 250 and t = 500. Figure 3.
Comparison of counter-propagating waves in the cRB equations (solid line) andthe constructed weakly-nonlinear solution (dashed line), at t = 150 (blue), t = 250 (red)and t = 500 (black), for (a) u and (b) w . Parameters are L = 300, N = 6000, ε = 0 . α = β = 1 . c = 1 . δ = γ = 1, k = 1 / √ k = (cid:112) . / x = 250, x = − t = 0 .
01 and ∆ T = ε ∆ t . weakly-nonlinear solutions. A further test was run by comparing the result of the interaction to acorresponding co-propagating case without interaction. These comparisons are shown in Figure 4for the left- and right-propagating radiating solitary waves independently. We can see that there isexcellent agreement between the cases, with only a very minor phase shift, suggesting that interactionof radiating solitary waves behaves essentially in the same way as solitary wave interaction. Now we consider the case when the characteristic speeds inthe equations are distinct, so we have c − O (1). For the case of co-propagating waves we takethe initial condition (4.1) and the results are presented in Figure 5. An Ostrovsky wave packet isformed in both layers, with the packet generated in the lower layer moving faster than the packet inthe upper layer. There is a reasonable agreement between the solutions, although the solution is lessaccurate away from the main wave packet due to the presence of radiation, which can re-enter thedomain. Furthermore, the solution for w is more accurate than the solution for u as this packet movesfaster and therefore will have some coupling in the upper layer for u , creating the disagreement. Thisis likely to be captured by the higher-order corrections, although the complexity of the calculationis increased by their inclusion. Alternatively, the errors can be reduced by decreasing the value of ε .The same case for ε = 0 . ε = 0 . n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -300 -200 -10000.511.5 -264 -26001 (a) Left-propagating solution for u at t =500.
100 200 30000.511.5 260 26401 (b) Right-propagating solution for u at t = 500. -300 -200 -10000.511.5 -264 -26000.51 (c) Left-propagating solution for w at t =500.
100 200 30000.511.5 260 26400.51 (d) Right-propagating solution for w at t = 500. Figure 4.
Comparison of the cRB equations in the counter-propagating case (blue, solid line)and the corresponding co-propagating case (red, dashed line), at t = 500 for (a, b) u , and (c,d) w . Parameters are L = 300, N = 6000, ε = 0 . α = β = 1 . c = 1 . δ = γ = 1, k = 1 / √ k = (cid:112) . / x = 250, x = − t = 0 .
01 and ∆ T = ε ∆ t . packet generated in the lower layer moving faster than the corresponding packet in the upper layer.At the point of interaction the packets move through each other and emerge with the appearanceof only small changes in their shape and structure, which could be attributed to evolution of thewave packet. This was confirmed by further simulations, which are excluded for brevity, showingthat interaction only causes a small phase shift in the main wave packet and some changes in thetail region due to fast-moving radiation. A second case of interest is using a cnoidal wave as the initial condition for the cRB equations. Wecan derive the cnoidal wave initial condition by considering the reduced case when u = w so that, in n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -200 0 2000123 -255 -250 -24502 (a) Solution for u at t = 0. -200 0 2000246 -255 -250 -24505 (b) Solution for w at t = 0. -200 0 200-10123 0 122.5 (c) Solution for u at t = 250. -200 0 200-20246 250 25155.5 (d) Solution for w at t = 250. Figure 5.
Comparison of co-propagating waves in the cRB equations (blue, solid line) andthe constructed weakly-nonlinear solution (red, dashed line), at t = 0, for (a) u and (b) w , and at t = 250 for (c) u and (d) w . Parameters are L = 300, N = 6000, ε = 0 . α = β = c = 2, δ = γ = 1, k = 1 / √ k = 1, x = 250, ∆ t = 0 .
01 and ∆ T = ε ∆ t . terms of w , we have w tt − c w xx = ε (cid:104) α (cid:0) w (cid:1) xx + βw ttxx (cid:105) . (4.3)The leading-order weakly-nonlinear solution to this equation takes the form of the KdV equation2 cf T + αf f ξ + βc f ξξξ = 0 , (4.4)where we have assumed that ξ = x − ct . This can be thought of as a reduction of the Ostrovsky casesderived in Section 2.2, under the assumption that the initial condition has zero mean. The cnoidalwave solution to this equation can be derived as f = − βc α (cid:32) f − ( f − f )cn (cid:34) ( ξ + νT ) (cid:114) f − f | m (cid:35)(cid:33) , (4.5) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -200 0 2000123 0 12.53 (a) Solution for u at t = 250. -200 0 2000246 250 2515.56 (b) Solution for w at t = 250. Figure 6.
Comparison of co-propagating waves in the cRB equations (blue, solid line) andthe constructed weakly-nonlinear solution (red, dashed line), at t = 250, for (a) u , and (b) w . Parameters are L = 300, N = 6000, ε = 0 . α = β = c = 2, δ = γ = 1, k = 1 / √ k = 1, x = 250, ∆ t = 0 .
01 and ∆ T = ε ∆ t . (a) Solution for u at multiple times. (b) Solution for w at multiple times. The wave movestwice as fast in this case. Figure 7.
Comparison of counter-propagating waves in the cRB equations (blue, solidline) and the constructed weakly-nonlinear solution (red, dashed line), for (a) u and (b) w . Parameters are L = 300, N = 6000, ε = 0 . α = β = c = 2, δ = γ = 1, k = 1 / √ k = 1, x = 250, x = − t = 0 .
01 and ∆ T = ε ∆ t . where ν = ( f + f + f ) βc, m = f − f f − f . n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction f < f < f and the restriction that 0 < m < L K = 2 K ( m ) (cid:114) f − f , (4.6)where K ( m ) is the complete elliptic integral of the first kind. This wave is an exact solution to theKdV equation and therefore will be an approximate solution to the cRB equations, to leading order,in the same way as the solitary wave solution in Section 4.1 was an approximate solution.To reduce the number of parameters in the solution, we will assume that f = 0, which meansthat the cnoidal wave is moving on a zero background. This will reduce the magnitude of the zeromean term, but it will not be zero, which is consistent with the approach for the solitary wave solution.To simplify notation we introduce the variable θ = (cid:113) f − f . Explicitly for the co-propagating wavewe have the initial condition F ( x ) = − f cn [( x + x ) θ | m ] , F ( x ) = − βc α f cn [( x + x ) θ | m ] ,V ( x ) = − f θ cn [( x + x ) θ | m ] sn [( x + x ) θ | m ] dn [( x + x ) θ | m ] ,V ( x ) = − βc α f θ cn [( x + x ) θ | m ] sn [( x + x ) θ | m ] dn [( x + x ) θ | m ] . (4.7) As was done for the solitary wave initial condition, let usconsider the case when the characteristic speeds in the equations are close, so we have c − O ( ε ).We take the initial condition for co-propagating waves, given in (4.7), and the results are presentedin Figure 8. Once again we see a good agreement between the results. As the troughs of the cnoidalwave are long, we can see the generation of a co-propagating radiating tail forming behind the cnoidalwave peaks, which then begins to interact with the preceding peak. The peaks of the wave appearto survive this interaction and maintain their shape. There is a reduction in amplitude but this canbe attributed to the generation of the radiating tail.The interacting case was also explored, with no visible change to the radiating solitary waves,as was seen for the solitary wave initial condition in the previous section. These are excluded as theimages appear qualitatively the same as the case of solitary wave initial conditions. Following on from Section 4.1, we now examine the case whenthe characteristic speeds in the equations are distinct, so we have c − O (1). We take the initialcondition (4.7) and the results are presented in Figure 9. As expected, an Ostrovsky wave packet isformed by each peak of the cnoidal wave, with the packets formed in the lower layer moving fasterthan the upper layer. The choice of initial condition here was such that the peaks of the cnoidal waveare distinct, so that each wave packet can be seen clearly. We can see that the agreement is goodaround the wave packet but becomes worse away from the main wave packet due to the radiationgenerated. This is consistent with the observation for the case of distinct characteristic speeds withsolitary wave initial conditions, however we note that the disagreement is much smaller in this case. n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -100 0 10000.20.40.60.81 94 98 102 00.51 (a) Solution for u at t = 0. -100 0 10000.20.40.60.81 94 98 102 00.51 (b) Solution for w at t = 0. -100 0 10000.20.40.60.81 106 110 00.51 (c) Solution for u at t = 300. -100 0 10000.20.40.60.81 106 110 0.51 (d) Solution for w at t = 300. Figure 8.
Comparison of co-propagating waves in the cRB equations (blue, solid line) andthe constructed weakly-nonlinear solution (red, dashed line), at t = 0, for (a) u and (b) w , and at t = 300 for (c) u and (d) w . Parameters are L ≈ . N = 2974, ε = 0 . f = 1 × − , f = 0, f = − , α = β = 1, c = 1 . δ = γ = 1, x = 0, ∆ t = 0 .
01 and∆ T = ε ∆ t . In the case of counter-propagating waves, due to the close proximity of the peaks for the cnoidalwave initial condition, the wave packets interact with each other and so we cannot observe the wavepacket evolution clearly. Therefore we exclude it from our considerations here.To see the evolution of the cnoidal wave initial condition into a series of Ostrovsky wave packets,we plot the evolution of the solution to the cRB equations over time in Figure 10, at time intervalsof t = 100 up to t = 1000. We have increased the value of δ and γ here to increase the speed of thepacket evolution. We can see that the peaks of the cnoidal wave gradually evolve into an Ostrovskywave packet and that the evolution occurs faster for u , as expected from the equations derived in n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -100 0 10000.20.40.60.81 94 98 102 00.51 (a) Solution for u at t = 0. -100 0 10001234 94 98 102 024 (b) Solution for w at t = 0. -100 0 100-0.500.5 100 105 -0.50 (c) Solution for u at t = 1000. -100 0 100-2-10123 108 110 112-1.6-1.4 (d) Solution for w at t = 1000. Figure 9.
Comparison of co-propagating waves in the cRB equations (blue, solid line) andthe constructed weakly-nonlinear solution (red, dashed line), at t = 0, for (a) u and (b) w ,and at t = 1000 for (c) u and (d) w . Parameters are L ≈ . N = 2974, ε = 0 . f = 1 × − , f = 0, f = − , α = β = 2, c = 2, δ = γ = 1, x = 0, ∆ t = 0 .
01 and∆ T = ε ∆ t . Section 2.2 due to the divisor of 2 c for w . Note that the domain is periodic, so although the waveappears to have a low speed, the time between images means that the previous peak is located closeto the current peak’s location. n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction -150 -100 -50 0 50 100 15002 0100 (a) Solution for u at various times. -150 -100 -50 0 50 100 15005 0100 (b) Solution for w at various times. Figure 10.
Solution to the cRB equations at intervals of t = 50, up to t = 1000, for (a) u , and (b) w . Parameters are L ≈ . N = 2974, ε = 0 . f = 1 × − , f = 0, f = − , α = β = 2, c = 2, δ = γ = 2, x = 0 and ∆ t = 0 . The system (1.1) - (1.2) is related to the system U tt − U xx = ε [ U x U xx + U ttxx − δ ( U − W )] , (5.1) W tt − c W xx = ε [ αW x W xx + βW ttxx + γ ( U − W )] , (5.2)by differentiation u = U x and w = W x . The system (5.1) - (5.2) is Lagrangian, and it has three localconservation laws for the mass, energy and momentum [10]: (cid:18) U t + δγ W t (cid:19) t − (cid:18) U x + δc γ W x + ε U x + εαδ γ W x + U ttx + εβδγ W ttx (cid:19) x = 0 , (5.3)12 (cid:26) U t + δγ W t + U x + δc γ W x + ε (cid:18) U x + αδγ W x (cid:19) + εU tx + εβδγ W tx + εδ ( U − W ) (cid:27) t − (cid:18) U t U x + δc γ W t W x + ε U t U x + εαδ γ W t W x + εU t U ttx + εβδγ W t W ttx (cid:19) x = 0 , (5.4) (cid:18) U t U x + δγ W t W x + εU tx U xx + εβδγ W tx W xx (cid:19) t − (cid:26) εU x U ttx + εβδγ W x W ttx + 12 (cid:20) U t + δγ W t + U x + δc γ W x + 2 ε (cid:18) U x + αδγ W x (cid:19) + εU tx + εβδγ W tx − εδ ( U − W ) (cid:21)(cid:27) x = 0 . (5.5) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction x and rewriting (5.4) and (5.5) in terms of u , w instead of U x , W x we obtain (cid:18) u t + δγ w t (cid:19) t − (cid:18) u x + δc γ w x + εuu x + εαδγ ww x + u ttx + εβδγ w ttx (cid:19) x = 0 , (5.6)12 (cid:26) U t + δγ W t + u + δc γ w + ε (cid:18) u + αδγ w (cid:19) + εu t + εβδγ w t + εδ ( U − W ) (cid:27) t − (cid:18) U t u + δc γ W t w + ε U t u + εαδ γ W t w + εU t u tt + εβδγ W t w tt (cid:19) x = 0 , (5.7) (cid:18) U t u + δγ W t w + εu t u x + εβδγ w t w x (cid:19) t − (cid:26) εuu tt + εβδγ ww tt + 12 (cid:20) U t + δγ W t + u + δc γ w + 2 ε (cid:18) u + αδγ w (cid:19) + εu t + εβδγ w t − εδ ( U − W ) (cid:21)(cid:27) x = 0 . (5.8)It is natural to assume that U ( − L ) = W ( − L ) = 0, and we make this assumption in what follows.Integrating these conservation laws with respect to x from − L to L , using the periodicity of u and w on [ − L, L ], we obtain one conserved quantity (mass)dd t (cid:18) ¯ u + δγ ¯ w (cid:19) = const , (5.9)and two non-local conservation laws (energy and momentum, respectively)12 dd t (cid:90) L − L (cid:26) U t + δγ W t + u + δc γ w + ε (cid:18) u + αδγ w (cid:19) + εu t + εβδγ w t + εδ ( U − W ) (cid:27) d x = (cid:26) U t ( L ) (cid:104) u ( L ) + ε u ( L ) + εu tt ( L ) (cid:105) + W t ( L ) (cid:20) δc γ w ( L ) + εαδ γ w ( L ) + εβδγ w tt ( L ) (cid:21)(cid:27) , (5.10)dd t (cid:90) L − L (cid:18) U t u + δγ W t w + εu t u x + εβδγ w t w x (cid:19) d x = 12 (cid:20) U t ( L ) + δγ W t ( L ) − εδ ( U ( L ) − W ( L )) (cid:21) , (5.11)where U t = ˆ U t + d¯ u d t , W t = ˆ W t + d ¯ w d t , (5.12)¯ u = 12 L (cid:90) L − L u d x = d + δd cos ωt, ¯ w = 12 L (cid:90) L − L w d x = d − γd cos ωt, (5.13)are the known periodic functions of time and ˆ U t , ˆ W t are zero-mean functions. As we impose that U ( − L ) = W ( − L ) = 0 it follows that U t ( − L ) = W t ( − L ) = 0. Therefore, we adjust the value of U t and W t to satisfy the zero boundary condition at x = − L .Note that if ¯ u = ¯ w = 0 (i.e. d = d = 0), then the nonlocal conservation laws (5.10) and(5.11) yield the usual (local) conservation of energy and momentum. More generally, the nonlocal n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction u and w .We now calculate the energy to confirm the derived relations (5.6) - (5.8). The parameters forthe simulation are chosen as ε = 0 . c = 1 .
025 and α = β = δ = γ = 1, with solitary wave initialconditions as taken in (3.1) with p = 1 to provide a significant non-zero mean value. The energyand momentum relations are plotted in Figure 11, where we note that the left-hand side of the massequation (5.9) is equal to zero in this case and therefore is omitted from the plot. (a) Energy for cRB equations. (b) Momentum for cRB equations. Figure 11.
Conserved quantities in the cRB equations, where the blue, solid line is the left-hand side of the relation and the red, dashed line is the right-hand side of the relation, for(a) energy and (b) momentum. Parameters are L = 300, N = 60000, c = 1 . ε = 0 . α = β = δ = γ = 1, k = k = 1 / √
2, ∆ t = 0 .
001 and ∆ T = ε ∆ t . We see that the energy and momentum oscillate with period f = 2 π/ω , although the energyrelation alternates between two different peaks. To verify the conservation laws, the peaks andtroughs were tracked for both the left-hand side and right-hand side of the conservation laws and theabsolute percentage error between peaks or troughs was calculated as 4 . × − % for the energyand 3 . × − % for the momentum. It is of note that the left-hand side of these laws require atime derivative on discrete data and so the accuracy to which they are conserved is hampered bythis limitation. In our case the time step was taken as ∆ t = 0 . In this paper we developed an asymptotic procedure for the construction of the d’Alembert-typesolution of the Cauchy problem for a system of coupled Boussinesq equations, which could describelong longitudinal strain waves in a bi-layer with an imperfect interface [10]. Two cases are presented,when the characteristic speeds in the equations are close or significantly different. An asymptoticexpansion is constructed to the case of significantly different characteristic speeds and the asymptoticexpansions were compared.We examined the accuracy of the constructed solution numerically, using direct numericalsimulations for the coupled Boussinesq equations and our constructed semi-analytical solution, and n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction u and w .The constructed solution can find useful applications in the studies of the scattering of radiatingsolitary waves by delamination [14], where the existing solutions do not take account of the zero-meanconstraint. The propagation of cnoidal waves in these structures is also of interest. Acknowledgements
MRT is grateful to the UK QJMAM Fund for Applied Mathematics for the support of his travel tothe ICoNSoM 2019 conference in Rome, Italy, where discussions of this work have taken place.
Appendix A Numerical Methods
To solve the equations derived in Section 2 and its subsections, we will make use of the pseudo-spectral methods for coupled Boussinesq equations in [25, 26] and single Ostrovsky equation in [21].For the coupled Ostrovsky equations we use the method outlined below, which also includes a briefoverview of the pseudo-spectral method.In the following method we use the Discrete Fourier Transform (DFT) to calculate the Fouriertransform of numerical data [27]. Let us consider a function u ( x, t ) on a finite domain x ∈ [ − L, L ]and we discretise the domain into N equally spaced points, so we have the spacing ∆ x = 2 L/N . Wescale the domain from x ∈ [ − L, L ] to ˜ x ∈ [0 , π ] via the transform ˜ x = sx + π , where s = π/L . n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction x j = − L + j ∆ x for j = 0 , . . . , N , we define the DFT for the function u ( x, t ) asˆ u ( k, t ) = 1 √ N N (cid:88) j =1 u ( x j , t ) e − ikx j , − N ≤ k ≤ N − , (A.1)and similarly the IDFT is defined as u ( x, t ) = 1 √ N N/ − (cid:88) k = − N/ ˆ u ( k, t ) e ikx j , j = 1 , , . . . , N, (A.2)where we have discretised and scaled wavenumber k ∈ Z . The transforms are implemented using theFFTW3 algorithm in C [28].We now consider the solution to the coupled Ostrovsky equations, where we can extend theRunge-Kutta method that was used for a single Ostrovsky equation in [21]. We present an examplefor the equation governing φ , which can be reduced to the equation for f . Explicitly we considerthe equations (2.29), (2.30), as this method can be reduced to solve the system (2.22). We changethe variables to x and t for simplicity, and introducing coefficients µ , ω for the terms in φ x and φ x respectively, we obtain(2 φ t + µφ x + ( f φ ) x + φ xxx ) x = δ ( φ − φ ) + H ( f ( x ) , f ( x )) , (2 φ t + ωφ x + α ( f φ ) x + βφ xxx ) x = γ ( φ − φ ) + H ( f ( x ) , f ( x )) , (A.3)where α , β , µ , ω , δ and γ are constants, and the functions f , f are known. We consider the equationon the domains t ∈ [0 , T ] and x ∈ [ − L, L ]. Here H , are functions of f , and are therefore known atthe given time step. Their form can be found in Section 2.1. Taking the Fourier transform of (A.3)gives (using the transform to ˜ x )2 ˆ φ t + (cid:0) iskµ − is k (cid:1) ˆ φ + isk F { f φ } = − iδsk (cid:16) ˆ φ − ˆ φ (cid:17) − isk ˆ H , φ t + (cid:0) iskω − is k β (cid:1) ˆ φ + iskα F { f φ } = − iγsk (cid:16) ˆ φ − ˆ φ (cid:17) − isk ˆ H . (A.4)Following the method of [27, 21] to remove the stiff term from this equation, we multiply through bya multiplicative factor M , and introduce a new function Φ , , where M , and Φ , take the form M = e − i ( s k − µsk − δsk ) t , ˆΦ = e − i ( s k − µsk − γsk ) t ˆ φ = M ˆ φ ,M = e − i ( βs k − ωsk − γsk ) t , ˆΦ = e − i ( βs k − ωsk − γsk ) t ˆ φ = M ˆ φ . (A.5)Substituting this into (A.4) leads to two ODEs for Φ , of the formˆΦ t = − isk M F (cid:40) F − (cid:34) ˆΦ M (cid:35)(cid:41) − i sk M (cid:18) ˆ S − δ Φ M (cid:19) , ˆΦ t = − iαsk M F (cid:40) F − (cid:34) ˆΦ M (cid:35)(cid:41) − i sk M (cid:18) ˆ S − γ Φ M (cid:19) . (A.6) n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction th order Runge-Kutta algorithm. We discretise the time domain as t j = j ∆ t and the functions as ˆΦ i,j = ˆΦ i ( k, t j ), ˆ φ i,j = ˆ φ i ( k, t j ) and ˆ f i,j = ˆ f i ( k, t j ), i = 1 ,
2. Introducing theadditional functions E = e i ( s k − µsk − δsk ) ∆ t , E = e i ( βs k − ωsk − γsk ) ∆ t , (A.7)we can introduce the optimised Runge-Kutta algorithm in the original variables φ , , taking the formˆ φ ,j +1 = E ˆ φ + 16 (cid:2) E k + 2 E ( k + k ) + k (cid:3) , ˆ φ ,j +1 = E ˆ φ + 16 (cid:2) E l + 2 E ( l + l ) + l (cid:3) ,k = − isk t F (cid:110) ˆ f ,j F − (cid:104) ˆ φ ,j (cid:105)(cid:111) − i ∆ t sk (cid:16) ˆ H − δ ˆ φ ,j (cid:17) ,l = − iαsk t F (cid:110) ˆ f ,j F − (cid:104) ˆ φ ,j (cid:105)(cid:111) − i ∆ t sk (cid:16) ˆ H − γ ˆ φ ,j (cid:17) ,k = − isk t F (cid:26) ˆ f ,j F − (cid:20) E (cid:18) ˆ φ ,j + k (cid:19)(cid:21)(cid:27) − i ∆ t sk (cid:18) ˆ H − δE (cid:20) φ ,j + l (cid:21)(cid:19) ,l = − iαsk t F (cid:26) ˆ f ,j F − (cid:20) E (cid:18) ˆ φ ,j + l (cid:19)(cid:21)(cid:27) − i ∆ t sk (cid:18) ˆ H − γE (cid:20) φ ,j + k (cid:21)(cid:19) ,k = − isk t F (cid:26) ˆ f ,j F − (cid:20) E ˆ φ ,j + k (cid:21)(cid:27) − i ∆ t sk (cid:18) ˆ H − δ (cid:20) E ˆ φ ,j + l (cid:21)(cid:19) ,l = − iαsk t F (cid:26) ˆ f ,j F − (cid:20) E ˆ φ ,j + l (cid:21)(cid:27) − i ∆ t sk (cid:18) ˆ H − γ (cid:20) E ˆ φ ,j + k (cid:21)(cid:19) ,k = − isk t F (cid:110) ˆ f ,j F − (cid:104) E ˆ φ ,j + E k (cid:105)(cid:111) − i ∆ t sk (cid:16) ˆ H − δ (cid:104) E ˆ φ ,j + E l (cid:105)(cid:17) ,l = − iαsk t F (cid:110) ˆ f ,j F − (cid:104) E ˆ φ ,j + E l (cid:105)(cid:111) − i ∆ t sk (cid:16) ˆ H − γ (cid:104) E ˆ φ ,j + E k (cid:105)(cid:17) . (A.8)When calculating the solution using this algorithm, the functions k i , l i must be calculated “in pairs”as the functions k and l are required when evaluating k and l , and so on. We can apply thisalgorithm to the case of a homogeneous Ostrovsky equation by setting ˆ H , = 0 and replacing theterm f i φ i with f i /
2, and therefore it can be applied to all derived coupled Ostrovsky equations inSection 2.2.We note that this is an extension of the method for a single Ostrovsky equation, as discussed in[21]. Therefore we can reduce the method for the coupled Ostrovsky equations to a single Ostrovskyequation as required.
References [1] G. A. Nariboli and A. Sedov. Burgers’s-Korteweg-De Vries equation for viscoelastic rods and plates.
J. Math.Anal. Appl. , 32:661–677, 1970.[2] L. A. Ostrovsky and A. M. Sutin. Nonlinear elastic waves in rods.
PMM , 41:531–537, 1977.[3] A. M. Samsonov.
Strain Solitons in Solids and How to Construct Them . Chapman & Hall/CRC, Boca Raton,2001.[4] A. V. Porubov.
Amplification of Nonlinear Strain Waves in Solids . World Scientific, Singapore, 2003. n Coupled Boussinesq Equations and Ostrovsky equations free from zero-mass contradiction [5] V. I. Erofeev, V. V. Kazhaev, and N. P. Semerikova. Waves in rods: dispersion, dissipation, nonlinearity .Fizmatlit, Moscow, 2002.[6] H.-H. Dai and X. Fan. Asymptotically approximate model equations for weakly nonlinear long waves incompressible elastic rods and their comparisons with other simplified model equations.
Math. Mech. Solids ,9:61–79, 2004.[7] T. Peets, K. Tamm, and J. Engelbrecht. On the role of nonlinearities in the Boussinesq-type wave equations.
Wave Motion , 71:113–119, 2017.[8] F. E. Garbuzov, K. R. Khusnutdinova, and I. V. Semenova. On Boussinesq-type models for long longitudinalwaves in elastic rods.
Wave Motion , 88:129–143, 2019.[9] F. E. Garbuzov, Y. M. Beltukov, and K. R. Khusnutdinova. Longitudinal bulk strain solitons in a hyperelasticrod with quadratic and cubic nonlinearities.
Theor. Math. Phys. , 202:319–333, 2020.[10] K. R. Khusnutdinova, A. M. Samsonov, and A. S. Zakharov. Nonlinear layered lattice model and generalizedsolitary waves in imperfectly bonded structures.
Phys. Rev. E , 79:056606, 2009.[11] K. R. Khusnutdinova and A. M. Samsonov. Fission of a longitudinal strain solitary wave in a delaminated bar.
Phys. Rev. E , 77:066603, 2008.[12] G. V. Dreiden, K. R. Khusnutdinova, A. M. Samsonov, and I. V. Semenova. Bulk strain solitary waves in bondedlayered polymeric bars with delamination.
J. Appl. Phys. , 112:063516, 2012.[13] K. R. Khusnutdinova and M. R. Tranter. Modelling of nonlinear wave scattering in a delaminated elastic bar.
Proc. Roy. Soc. A , 471:20150584, 2015.[14] K. R. Khusnutdinova and M. R. Tranter. On radiating solitary waves in bi-layers with delamination and coupledOstrovsky equations.
Chaos , 27:013112, 2017.[15] A. V. Belashov, Y. M. Beltukov, and I. V. Semenova. Pump-probe digital holography for monitoring of long bulknonlinear strain waves in solid waveguides.
Proc. SPIE , 10678:1067810, 2018.[16] N. Peake and S. V. Sorokin. A nonlinear model of the dynamics of a large elastic plate with heavy fluid loading.
Proc. R. Soc. A , 462:2205–2224, 2006.[17] D. A. Indejtsev, M. G. Zhuchkova, D. P. Kouzov, and S. V. Sorokin. Low-frequency wave propagation in anelastic plate floating on a two-layered fluid.
Wave Motion , 62:98–113, 2016.[18] Z. Abiza, M. Destrade, and R. W. Ogden. Large acoustoelastic effect.
Wave Motion , 49:364–374, 2012.[19] I. V. Andrianov, V. D. Danishevsky, J. D. Kaplunov, and B. Markert. Wide frequency higher-order dynamicmodel for transient waves in a lattice. In I. V. Andrianov, A. I. Manevich, Y. V. Mikhlin, and O. V. Gendelman,editors,
Problems of Nonlinear Mechanics and Physics of Materials , chapter 1, pages 3–12. Springer, 2019.[20] K. R. Khusnutdinova, K. R. Moore, and D. E. Pelinovsky. Validity of the weakly nonlinear solution of the Cauchyproblem for the Boussinesq-type equation.
Stud. Appl. Math. , 133:52–83, 2014.[21] K. R. Khusnutdinova and M. R. Tranter. D’Alembert-type solution of the Cauchy problem for the Boussinesq-Klein-Gordon equation.
Stud. Appl. Math. , 142:551–585, 2019.[22] L. A. Ostrovsky. Nonlinear internal waves in a rotating ocean.
Oceanology , 18:119–125, 1978.[23] K. R. Khusnutdinova and K. R. Moore. Initial-value problem for coupled Boussinesq equations and a hierarchyof Ostrovsky equations.
Wave Motion , 48:738–752, 2011.[24] R. H. J. Grimshaw, K. R. Khusnutdinova, and K. R. Moore. Radiating solitary waves in coupled Boussinesqequations.
IMA J. Appl. Math. , 82:802–820, 2017.[25] K. R. Khusnutdinova and M. R. Tranter. Weakly-nonlinear solution of coupled Boussinesq equations and radiatingsolitary waves. In H. Altenbach, A. Belyaev, V. Eremeyev, A. Krivtsov, and A. Porubov, editors,
DynamicalProcesses in Generalized Continua and Structures , chapter 18, pages 321–343. Springer, 2019.[26] J. Engelbrecht, A. Salupere, and K. Tamm. Waves in microstructured solids and the Boussinesq paradigm.
WaveMotion , 48:717–726, 2011.[27] L. N. Trefethen.
Spectral Methods in MATLAB . SIAM, Philadelphia, 2000.[28] M. Frigo and S. G. Johnson. The design and implementation of FFTW3.