# Numerical Solution of the Modified Equal Width Wave Equation

Authors: Seydi Battal Gazi Karakoç; Turabi Geyikli
Journal: International Journal of Differential Equations
 (2012)

---

## Abstract

Numerical solution of the modified equal width wave equation is obtained by using lumped Galerkin method based on cubic B-spline finite element method. Solitary wave motion and interaction of two solitary waves are studied using the proposed method. Accuracy of the proposed method is discussed by computing the numerical conserved lawsL2 and L∞ error norms. The numerical results are found in good agreement with exact solution. A linear stability analysis of the scheme is also investigated.


---

## Body

## 1. Introduction

The modified equal width wave equation (MEW) based upon the equal width wave (EW) equation [1, 2] which was suggested by Morrison et al. [3] is used as a model partial differential equation for the simulation of one-dimensional wave propagation in nonlinear media with dispersion processes. This equation is related with the modified regularized long wave (MRLW) equation [4] and modified Korteweg-de Vries (MKdV) equation [5]. All the modified equations are nonlinear wave equations with cubic nonlinearities and all of them have solitary wave solutions, which are wave packets or pulses. These waves propagate in non-linear media by keeping wave forms and velocity even after interaction occurs. Few analytical solutions of the MEW equation are known. Thus numerical solutions of the MEW equation can be important and comparison between analytic solution can be made. Geyikli and Battal Gazi Karakoç [6, 7] solved the MEW equation by a collocation method using septic B-spline finite elements and using a Petrov-Galerkin finite element method with weight functions quadratic and element shape functions which are cubic B-splines. Esen applied a lumped Galerkin method based on quadratic B-spline finite elements which have been used for solving the EW and MEW equation [8, 9]. Saka proposed algorithms for the numerical solution of the MEW equation using quintic B-spline collocation method [10]. Zaki considered the solitary wave interactions for the MEW equation by collocation method using quintic B-spline finite elements [11] and obtained the numerical solution of the EW equation by using least-squares method [12]. Wazwaz investigated the MEW equation and two of its variants by the tanh and the sine-cosine methods [13]. A solution based on a collocation method incorporated cubic B-splines is investigated by and Saka and Dağ [14]. Variational iteration method is introduced to solve the MEW equation by Lu [15]. Evans and Raslan [16] studied the generalized EW equation by using collocation method based on quadratic B-splines to obtain the numerical solutions of a single solitary waves and the birth of solitons. Hamdi et al. [17] derived exact solitary wave solutions of the generalized EW equation using Maple software. Esen and Kutluay studied a linearized implicit finite difference method in solving the MEW equation [18]. In the present work we solve the MEW equation numerically by a lumped Galerkin method using cubic B-spline finite elements. The accuracy of the proposed method is demonstrated by two test problems: the motion of a single solitary wave and the interaction of two solitary waves. A linear stability analysis based on a Fourier method shows that the numerical scheme is unconditionally stable.


## 2. Cubic B-Spline Lumped Galerkin Method

The modified equal width wave (MEW) equation considered here has the normalized form [3](2.1)Ut+3U2Ux-μUxxt=0,
with the physical boundary conditions U→0 as x→±∞, where t is time, x is the space coordinate, μ is a positive parameter, and U(x,t) is wave amplitude. In this study, boundary conditions are chosen from(2.2)U(a,t)=0,U(b,t)=0,Ux(a,t)=0,Ux(b,t)=0,t>0,
and the initial condition(2.3)U(x,0)=f(x),a≤x≤b,
where f(x) is a localized disturbance inside interval [a,b]. The interval [a,b] is partitioned into uniformly sized finite elements of length h by the knots xm such that a=x0<x1⋯<xN=b and h=(xm+1-xm). The cubic B-splines ϕm(x), (m=-1(1)N+1), at the knots xm are defined over the interval [a,b] by [19](2.4)ϕm(x)=1h3{(x-xm-2)3,x∈[xm-2,xm-1],h3+3h2(x-xm-1)+3h(x-xm-1)2-3(x-xm-1)3,x∈[xm-1,xm],h3+3h2(xm+1-x)+3h(xm+1-x)2-3(xm+1-x)3,x∈[xm,xm+1],(xm+2-x)3,x∈[xm+1,xm+2],0,otherwise.
The set of functions {ϕ-1(x),ϕ0(x),…,ϕN+1(x)} forms a basis for functions defined over [a,b]. The approximate solution UN(x,t) to the exact solution U(x,t) is given by(2.5)UN(x,t)=∑j=-1N+1ϕj(x)δj(t),
where δj are time-dependent parameters to be determined from the boundary and weighted residual conditions. Each cubic B-spline covers 4 elements so that each element [xm,xm+1] is covered by 4 splines. In each element, using the following local coordinate transformation(2.6)hξ=x-xm,0≤ξ≤1,
cubic B-spline shape functions in terms of ξ over the element [0,1] can be defined as(2.7)ϕm-1ϕmϕm+1ϕm+2={(1-ξ)3,1+3(1-ξ)+3(1-ξ)2-3(1-ξ)3,1+3ξ+3ξ2-3ξ3,ξ3.
All splines apart from ϕm-1(x),ϕm(x),ϕm+1(x), and ϕm+2(x) are zero over the element [0,1]. Variation of the function U(x,t) over element [0,1] is approximated by(2.8)UN(ξ,t)=∑j=m-1m+2δjϕj,
where δm-1,δm,δm+1,δm+2 act as element parameters and B-splines ϕm-1,ϕm,ϕm+1,ϕm+2 as element shape functions. Using trial function (2.5) and cubic splines (2.4), the values of U,U′,U′′ at the knots are determined in terms of the element parameters δm by(2.9)Um=U(xm)=δm-1+4δm+δm+1,Um′=U′(xm)=3(-δm-1+δm+1),Um′′=U′′(xm)=6(δm-1-2δm+δm+1),
where the symbols ′ and ′′ denote first and second differentiation with respect to x, respectively. The splines ϕm(x) and its two principle derivatives vanish outside the interval [xm-2,xm+2]. Use Galerkin’s method with weight function W(x) to obtain the weak form of (2.1) which is(2.10)∫abW(Ut+3U2Ux-μUxxt)dx=0.
For a single element [xm,xm+1] using transformation (2.6) into the (2.10) we obtain
(2.11)∫01W(Ut+3hÛ2Uξ-μh2Uξξt)dξ=0,
where Û is taken to be a constant over an element to simplify the integral. Integrating (2.11) by parts leads to (2.12)∫01[WUt+λWUξ+βWξUξt]dξ=βWUξt|01
where λ=3Û2/h and β=μ/h2. Taking the weight function with cubic B-spline shape functions given by (2.7) and substituting approximation (2.8) in integral equation (2.12) with some manipulation, we obtain the element contributions in the form(2.13)∑j=m-1m+2[(∫01ϕiϕj+βϕi′ϕj′)dξ-βϕiϕj′|01]δ̇je+∑j=m-1m+2(λ∫01ϕiϕj′dξ)δje.
In matrix notation this equation becomes(2.14)[Ae+β(Be-Ce)]δ̇e+λDeδe,
where δe=(δm-1,δm,δm+1,δm+2)T are the element parameters and the dot denotes differentiation with respect to t. The element matrices Ae,Be,Ce, and De are given by the following integrals:(2.15)Aije=∫01ϕiϕjdξ=1140[2012960112911889336060933118812916012920],Bije=∫01ϕi′ϕj′dξ=110[1821-36-321102-87-36-36-8710212-3-362118],Cije=ϕiϕj′|01=3[10-104-1-411-4-140-101],Dije=∫01ϕiϕj′dξ=120[-10-9181-71-15018338-38-18315071-1-18910],
where the suffices i,j take only the values m-1,m,m+1,m+2 for the typical element [xm,xm+1]. A lumped value for λ is found from (1/4)(Um+Um+1)2 as(2.16)λ=34h(δm-1+5δm+5δm+1+δm+2)2.
By assembling all contributions from all elements, (2.14) leads to the following matrix equation:(2.17)[Ae+β(Be-Ce)]δ̇e+λDeδe=0,
where δ=(δ-1,δ0⋯δN,δN+1)T is a global element parameter. The matrices A,B, and λD are septadiagonal and row of each has the following form:(2.18)A=1140(1,120,1191,2416,1191,120,1),B=110(-3,-72,-45,240,-45,-72,-3),λD=120(-λ1,-18λ1-38λ2,9λ1-183λ2-71λ3,10λ1+150λ2-150λ3-10λ4,71λ2+183λ3-9λ4,38λ3+18λ4,λ4),
where(2.19)λ1=34h(δm-2+5δm-1+5δm+δm+1)2,λ2=34h(δm-1+5δm+5δm+1+δm+2)2,λ3=34h(δm+5δm+1+5δm+2+δm+3)2,λ4=34h(δm+1+5δm+2+5δm+3+δm+4)2.
Replacing the time derivative of the parameter δ̇ by usual forward finite difference approximation and parameter δ by the Crank-Nicolson formulation(2.20)δ̇=δn+1-δnΔt,δ=12(δn+δn+1)
into equation (2.17), gives the (N+3)×(N+3) septadiagonal matrix system(2.21)[A+β(B-C)+λΔt2D]δn+1=[A+β(B-C)-λΔt2D]δn,
where Δt is time step. Applying the boundary conditions (2.2) to the system (2.21) we obtain an (N+1)×(N+1) septadiagonal matrix system. This system is efficiently solved with a variant of the Thomas algorithm, but an inner iteration is also needed at each time step to cope with the non-linear term. A typical member of the matrix system (2.21) may be written in terms of the nodal parameters δn and δn+1 as(2.22)γ1δm-2n+1+γ2δm-1n+1+γ3δmn+1+γ4δm+1n+1+γ5δm+2n+1+γ6δm+3n+1+γ7δm+4n+1=γ7δm-2n+γ6δm-1n+γ5δmn+γ4δm+1n+γ3δm+2n+γ2δm+3n+γ1δm+4n,
where(2.23)γ1=1140-3β10-λΔt40,γ2=120140-72β10-56λΔt40,γ3=1191140-45β10-245λΔt40,γ4=2416140+240β10,γ5=1191140-45β10+245λΔt40,γ6=120140-72β10+56λΔt40,γ7=1140-3β10+λΔt40
which all depend on δn. The initial vector of parameter δ0=(δ-10,…,δN+10) must be determined to iterate system (2.21). To do this, the approximation is rewritten over the interval [a,b] at time t=0 as follows:(2.24)UN(x,0)=∑m=-1N+1ϕm(x)δm0,
where the parameters δm0 will be determined. UN(x,0) are required to satisfy the following relations at the mesh points xm:(2.25)UN(xm,0)=U(xm,0),m=0,1,…,N,UN′(x0,0)=U′(xN,0)=0.
The above conditions lead to a tridiagonal matrix system of the form(2.26)[-303141⋱141-303][δ-01δ00⋮δN0δN+10]=[0U(x0)⋮U(xN)0]
which can be solved using a variant of the Thomas algorithm.


### 2.1. Stability Analysis

The stability analysis is based on the Von Neumann theory in which the growth factorg of the error in a typical mode of amplitude δ̂n,(2.27)δjn=δ̂neijkh,
where k is the mode number and h the element size, is determined from a linearization of the numerical scheme. To apply the stability analysis, the MEW equation needs to be linearized by assuming that the quantity U in the non-linear term U2Ux is locally constant. Substituting the Fourier mode (2.27) into (2.22) gives the growth factor g of the form(2.28)g=a-iba+ib,
where(2.29)a=2416+3360β+(2382-1260β)cos⁡θh+(240-2016β)cos⁡2θh+(2-84β)cos⁡3θh,b=5145λΔtsin⁡θh+1176λΔtsin⁡2θh+21λΔtsin⁡3θh.
The modulus of |g| is 1, therefore the linearized scheme is unconditionally stable.


## 2.1. Stability Analysis

The stability analysis is based on the Von Neumann theory in which the growth factorg of the error in a typical mode of amplitude δ̂n,(2.27)δjn=δ̂neijkh,
where k is the mode number and h the element size, is determined from a linearization of the numerical scheme. To apply the stability analysis, the MEW equation needs to be linearized by assuming that the quantity U in the non-linear term U2Ux is locally constant. Substituting the Fourier mode (2.27) into (2.22) gives the growth factor g of the form(2.28)g=a-iba+ib,
where(2.29)a=2416+3360β+(2382-1260β)cos⁡θh+(240-2016β)cos⁡2θh+(2-84β)cos⁡3θh,b=5145λΔtsin⁡θh+1176λΔtsin⁡2θh+21λΔtsin⁡3θh.
The modulus of |g| is 1, therefore the linearized scheme is unconditionally stable.


## 3. Numerical Examples and Results

In this part, we consider the following two test problems: the motion of a single solitary wave and interaction of two solitary waves. All computations are executed on a Pentium 4 PC in the Fortran code using double precision arithmetic. For the MEW equation, it is important to discuss the following three invariant conditions given in [11], which, respectively, correspond to conversation of mass, momentum, and energy. The accuracy of the method is measured by both the L2 error norm(3.1)C1=∫abUdx≃h∑J=1NUjn,C2=∫abU2+μ(Ux)2dx≃h∑J=1N(Ujn)2+μ(Ux)jn,C3=∫abU4dx≃h∑J=1N(Ujn)4,(3.2)L2=‖Uexact-UN‖2≃h∑J=0N|Ujexact-(UN)j|2
and the L∞ error norm(3.3)L∞=‖Uexact-UN‖∞≃max⁡j|Ujexact-(UN)j|
to show how well the scheme predicts the position and amplitude of the solution as the simulation proceeds. The variable Ujn and its first derivative appearing in (3.1) can be computed from (2.9), respectively.


### 3.1. The Motion of Single Solitary Wave

For this problem we consider (2.1) with the boundary condition U→0 as x→±∞ and the initial condition(3.4)U(x,0)=Asech[k(x-x0)].
An exact solution of this problem is given by [11](3.5)U(x,t)=Asech[k(x-x0-vt)]
which represents the motion of a single solitary wave with amplitude A, where the wave velocity v=A2/2 and k=1/μ. For this problem the analytical values of the invariants are [11](3.6)C1=Aπk,C2=2A2k+2μkA23,C3=4A43k.
For the numerical simulation of the motion of a single solitary wave, we choose the parameters h=0.1, Δt=0.05, μ=1, x0=30, A=0.25 through the interval 0≤x≤80. The analytical values for the invariants are C1=0.7853982, C2=0.1666667, C3=0.0052083. The invariants C1 and C2 change from their initial values by less than 1×10-7 and 2×10-7 respectively, during the time of running, whereas the changes of invariant C3 approach to zero throughout. The computations are done until time t=20, and we find L2, L∞ error norms and numerical invariants C1, C2, C3 at various times. Results are documented in Table 1. One may also compare our results with those in the other studies [9, 16, 18]. According to both L2, L∞ error norms, agreement between numerical values and exact solution appears very satisfactorily through illustrations of three invariants and norms. Figure 1 shows that the proposed method performs the motion of propagation of a solitary wave satisfactorily, which moved to the right at a constant speed and preserved its amplitude and shape with increasing time as expected. Amplitude is 0.25 at t=0 which is located at x=30, while it is 0.249900 at t=20 which is located at x=30.6. The absolute difference in amplitudes at times t=0 and t=20 is 1×10-4 so that there is a little change between amplitudes. The convergence rates for the numerical method in space sizes hm and time steps UΔtm can be calculated by following formula [9], respectively,(3.7)order=log10(|Uexact-Uhmnum|/|Uexact-Uhm+1num|)log10(hm/hm+1),order=log10(|Uexact-UΔtmnum|/|Uexact-UΔtm+1num|)log10(Δtm/Δtm+1).
Table 2 displays the convergence rates for different values of space size h and a fixed value of the time step Δtm. We have clearly seen that the convergence rates when Δt is fixed are not good these for size. In addition, the time rate of the convergence for the numerical method is computed with various time step UΔtm and fixed space step h in Table 3. It can clearly be seen that the present method provides remarkable reductions in convergence rates for the smaller time.Table 1
Invariants and error norms for single solitary wave withh=0.1, Δt=0.05.
tC1C2C3L2×103L∞×10300.78539660.16666610.00520830.00000000.000000050.78539660.16666620.00520830.02048380.0115451100.78539660.16666620.00520830.04077430.0231561150.78539670.16666620.00520830.06069750.0347169200.78539670.16666630.00520830.08009800.046061820 [9]0.78538980.16676140.00520820.07969400.046552320 [16]0.78495450.16647650.00519950.29051660.249892520 [18]0.78539770.16647350.00520830.26928120.2569972Table 2
The order of convergence att=20,Δt=0.05,A=0.25.
hmL2×103orderL∞×103order0.84.16296467—2.78665608—0.41.212289311.779877270.684622042.025155310.20.317526401.932785580.181756451.913301170.10.080098011.98038240.046061811.980363550.050.019324482.051336800.011229742.036246570.0250.005309591.863757220.003040441.88497250Table 3
The order of convergence att=20, h=0.1, A=0.25.
ΔtmL2×103orderL∞×103order0.80.08383192—0.08300304—0.40.074244890.175207930.050611520.713698370.20.07833790−0.078175430.044485030.186145870.10.07972878−0.025390130.04573242−0.039897330.050.08009801−0.006665800.04606181−0.010353830.0250.08019167−0.001685980.04614408−0.00257446The motion of a single solitary wave withh=1,Δt=0.05 at t=0 and t=20.
(a)(b)


### 3.2. Interaction of Two Solitary Waves

In this section, we consider (2.1) with boundary conditions U→0 as x→±∞ and the initial condition(3.8)U(x,0)=∑j=12Ajsech⁡(k[x-xj])
where k=1/μ.Firstly we studied the interaction of two positive solitary waves with the parametersh=0.1, Δt=0.025, μ=1, A1=1, A2=0.5, x1=15, x2=30 through the interval 0≤x≤80 which used the earlier papers [9–11]. The analytic invariants are [16] C1=π(A1+A2)=4.7123889, C2=(8/3)(A12+A22)=3.3333333, C3=(4/3)(A14+A24)=1.4166667 and changes in C1,C2, and C3 are less than 4.1×10-3, 4.3×10-3, and 3.6×10-3 percent, respectively, as can be seen in Table 4. The experiment was run from t=0 to t=55 to allow the interaction to take place. This condition is propagated to the right with velocities dependent upon their magnitudes and a stage is reached where the larger wave has passed through the smaller solitary wave and has emerged with their original shapes. Figure 2 shows the interactions of two positive solitary waves. Interaction started at about time t=25, overlapping processes occurred between times t=25 and t=40 and, waves started to resume their original shapes after time t=40. It can be seen that, at t=5, the wave with larger amplitude is on the left of the second wave with smaller amplitude. The larger wave catches up with the smaller one as time increases. At t=55, the amplitude of larger waves is 0.999581 at the point x=44.4 whereas the amplitude of the smaller one is 0.510464 at the point x=34.7. It is found that the absolute difference in amplitude is 1.04×10-1 for the smaller wave and 0.419×10-3 for the larger wave for this algorithm. Finally, we have studied the interaction of two solitary waves with the following parameters, μ=1, x1=15, x2=30, A1=-2, A2=1 together with time step Δt=0.025 and space step h=0.1 in the range 0≤x≤150. The experiment was run from t=0 to t=55 to allow the interaction to take place. Figure 3 shows the development of the solitary wave interaction. As is seen from Figure 3, at t=0, a wave with the negative amplitude is on the left of another wave with the positive amplitude. The larger wave with the negative amplitude catches up with the smaller one with the positive amplitude as the time increases. At t=55, the amplitude of the smaller wave is 0.972910 at the point x=52.5, whereas the amplitude of the larger one is −2.016990 at the point x=124.3. It is found that the absolute difference in amplitudes is 0.27×10-1 for the smaller wave and 0.16×10-1 for the larger wave. The analytical invariants can be found as C1=-3.1415927,C2=13.3333333,C3=22.6666667 and changes in C1,C2, and C3 are less than 10.5×10-3, 3.7×10-3, and 10.9×10-3, percent, respectively. Table 4 lists the values of the invariants of the two solitary waves with amplitudes A1=1,A2=0.5, and A1=-2,A2=1. It can be seen that the values obtained for the invariants are satisfactorily constant during the computer run. We have also compared the computed values of the invariants of the two solitary waves with results from [9] in Table 5.Table 4
Invariants for the interaction of two solitary waves.
A1=1,A2=0.5(0≤x≤80)A1=-2,A2=1(0≤x≤150)tC1C2C3C1C2C304.71237323.33332531.4166643−3.141573713.333281622.666531354.71238613.33334821.4166852−3.145860313.344984322.7133525104.71239593.33336211.4166982−3.137754313.303115322.5832812154.71240653.33337851.4167141−3.162543613.383899122.8926223204.71242493.33341641.4167521−3.165831813.399930422.9451481254.71248993.33358321.4169238−3.170181913.412639122.9954601304.71276433.33335571.4177617−3.174755313.425135822.0458834354.71304743.33525001.4188849−3.179370713.437678523.0966349404.71248813.33363161.4171690−3.184012613.450289523.1477374454.71230023.33318781.4167580−3.188678913.462973023.1991979504.71224793.33309231.4167142−3.193369913.475730323.2510209554.71225763.33311491.4167237−3.198085613.488562423.3032108Table 5
Comparison of invariants for the interaction of two solitary waves with results from [9] with h=0.1, Δt=0.025, A1=1, A2=0.5, (0≤x≤80).
tC1C2C3C1[9]C2[9]C3[9]04.71237323.33332531.41666434.71238843.33528901.416669754.71238613.33334821.41668524.71237183.33526351.4166486104.71239593.33336211.41669824.71238533.33528361.4166647154.71240653.33337851.41671414.71237563.33528941.4166772204.71242493.33341641.41675214.71237483.33530411.4166926254.71248993.33358321.41692384.71241733.33542781.4168363304.71276433.33335571.41776174.71264103.33594641.4176398354.71304743.33525001.41888494.71283533.33642471.4186746404.71248813.33363161.41716904.71239463.33559511.4170695454.71230023.33318781.41675804.71222733.33523641.4166637504.71224793.33309231.41671424.71215673.33511751.4165797554.71225763.33311491.41672374.71214003.33508471.4165527Interaction of two solitary waves at different times.
(a)(b)(c)(d)(e)(f)Interaction of two solitary waves at different times.
(a)(b)(c)(d)(e)(f)


## 3.1. The Motion of Single Solitary Wave

For this problem we consider (2.1) with the boundary condition U→0 as x→±∞ and the initial condition(3.4)U(x,0)=Asech[k(x-x0)].
An exact solution of this problem is given by [11](3.5)U(x,t)=Asech[k(x-x0-vt)]
which represents the motion of a single solitary wave with amplitude A, where the wave velocity v=A2/2 and k=1/μ. For this problem the analytical values of the invariants are [11](3.6)C1=Aπk,C2=2A2k+2μkA23,C3=4A43k.
For the numerical simulation of the motion of a single solitary wave, we choose the parameters h=0.1, Δt=0.05, μ=1, x0=30, A=0.25 through the interval 0≤x≤80. The analytical values for the invariants are C1=0.7853982, C2=0.1666667, C3=0.0052083. The invariants C1 and C2 change from their initial values by less than 1×10-7 and 2×10-7 respectively, during the time of running, whereas the changes of invariant C3 approach to zero throughout. The computations are done until time t=20, and we find L2, L∞ error norms and numerical invariants C1, C2, C3 at various times. Results are documented in Table 1. One may also compare our results with those in the other studies [9, 16, 18]. According to both L2, L∞ error norms, agreement between numerical values and exact solution appears very satisfactorily through illustrations of three invariants and norms. Figure 1 shows that the proposed method performs the motion of propagation of a solitary wave satisfactorily, which moved to the right at a constant speed and preserved its amplitude and shape with increasing time as expected. Amplitude is 0.25 at t=0 which is located at x=30, while it is 0.249900 at t=20 which is located at x=30.6. The absolute difference in amplitudes at times t=0 and t=20 is 1×10-4 so that there is a little change between amplitudes. The convergence rates for the numerical method in space sizes hm and time steps UΔtm can be calculated by following formula [9], respectively,(3.7)order=log10(|Uexact-Uhmnum|/|Uexact-Uhm+1num|)log10(hm/hm+1),order=log10(|Uexact-UΔtmnum|/|Uexact-UΔtm+1num|)log10(Δtm/Δtm+1).
Table 2 displays the convergence rates for different values of space size h and a fixed value of the time step Δtm. We have clearly seen that the convergence rates when Δt is fixed are not good these for size. In addition, the time rate of the convergence for the numerical method is computed with various time step UΔtm and fixed space step h in Table 3. It can clearly be seen that the present method provides remarkable reductions in convergence rates for the smaller time.Table 1
Invariants and error norms for single solitary wave withh=0.1, Δt=0.05.
tC1C2C3L2×103L∞×10300.78539660.16666610.00520830.00000000.000000050.78539660.16666620.00520830.02048380.0115451100.78539660.16666620.00520830.04077430.0231561150.78539670.16666620.00520830.06069750.0347169200.78539670.16666630.00520830.08009800.046061820 [9]0.78538980.16676140.00520820.07969400.046552320 [16]0.78495450.16647650.00519950.29051660.249892520 [18]0.78539770.16647350.00520830.26928120.2569972Table 2
The order of convergence att=20,Δt=0.05,A=0.25.
hmL2×103orderL∞×103order0.84.16296467—2.78665608—0.41.212289311.779877270.684622042.025155310.20.317526401.932785580.181756451.913301170.10.080098011.98038240.046061811.980363550.050.019324482.051336800.011229742.036246570.0250.005309591.863757220.003040441.88497250Table 3
The order of convergence att=20, h=0.1, A=0.25.
ΔtmL2×103orderL∞×103order0.80.08383192—0.08300304—0.40.074244890.175207930.050611520.713698370.20.07833790−0.078175430.044485030.186145870.10.07972878−0.025390130.04573242−0.039897330.050.08009801−0.006665800.04606181−0.010353830.0250.08019167−0.001685980.04614408−0.00257446The motion of a single solitary wave withh=1,Δt=0.05 at t=0 and t=20.
(a)(b)


## 3.2. Interaction of Two Solitary Waves

In this section, we consider (2.1) with boundary conditions U→0 as x→±∞ and the initial condition(3.8)U(x,0)=∑j=12Ajsech⁡(k[x-xj])
where k=1/μ.Firstly we studied the interaction of two positive solitary waves with the parametersh=0.1, Δt=0.025, μ=1, A1=1, A2=0.5, x1=15, x2=30 through the interval 0≤x≤80 which used the earlier papers [9–11]. The analytic invariants are [16] C1=π(A1+A2)=4.7123889, C2=(8/3)(A12+A22)=3.3333333, C3=(4/3)(A14+A24)=1.4166667 and changes in C1,C2, and C3 are less than 4.1×10-3, 4.3×10-3, and 3.6×10-3 percent, respectively, as can be seen in Table 4. The experiment was run from t=0 to t=55 to allow the interaction to take place. This condition is propagated to the right with velocities dependent upon their magnitudes and a stage is reached where the larger wave has passed through the smaller solitary wave and has emerged with their original shapes. Figure 2 shows the interactions of two positive solitary waves. Interaction started at about time t=25, overlapping processes occurred between times t=25 and t=40 and, waves started to resume their original shapes after time t=40. It can be seen that, at t=5, the wave with larger amplitude is on the left of the second wave with smaller amplitude. The larger wave catches up with the smaller one as time increases. At t=55, the amplitude of larger waves is 0.999581 at the point x=44.4 whereas the amplitude of the smaller one is 0.510464 at the point x=34.7. It is found that the absolute difference in amplitude is 1.04×10-1 for the smaller wave and 0.419×10-3 for the larger wave for this algorithm. Finally, we have studied the interaction of two solitary waves with the following parameters, μ=1, x1=15, x2=30, A1=-2, A2=1 together with time step Δt=0.025 and space step h=0.1 in the range 0≤x≤150. The experiment was run from t=0 to t=55 to allow the interaction to take place. Figure 3 shows the development of the solitary wave interaction. As is seen from Figure 3, at t=0, a wave with the negative amplitude is on the left of another wave with the positive amplitude. The larger wave with the negative amplitude catches up with the smaller one with the positive amplitude as the time increases. At t=55, the amplitude of the smaller wave is 0.972910 at the point x=52.5, whereas the amplitude of the larger one is −2.016990 at the point x=124.3. It is found that the absolute difference in amplitudes is 0.27×10-1 for the smaller wave and 0.16×10-1 for the larger wave. The analytical invariants can be found as C1=-3.1415927,C2=13.3333333,C3=22.6666667 and changes in C1,C2, and C3 are less than 10.5×10-3, 3.7×10-3, and 10.9×10-3, percent, respectively. Table 4 lists the values of the invariants of the two solitary waves with amplitudes A1=1,A2=0.5, and A1=-2,A2=1. It can be seen that the values obtained for the invariants are satisfactorily constant during the computer run. We have also compared the computed values of the invariants of the two solitary waves with results from [9] in Table 5.Table 4
Invariants for the interaction of two solitary waves.
A1=1,A2=0.5(0≤x≤80)A1=-2,A2=1(0≤x≤150)tC1C2C3C1C2C304.71237323.33332531.4166643−3.141573713.333281622.666531354.71238613.33334821.4166852−3.145860313.344984322.7133525104.71239593.33336211.4166982−3.137754313.303115322.5832812154.71240653.33337851.4167141−3.162543613.383899122.8926223204.71242493.33341641.4167521−3.165831813.399930422.9451481254.71248993.33358321.4169238−3.170181913.412639122.9954601304.71276433.33335571.4177617−3.174755313.425135822.0458834354.71304743.33525001.4188849−3.179370713.437678523.0966349404.71248813.33363161.4171690−3.184012613.450289523.1477374454.71230023.33318781.4167580−3.188678913.462973023.1991979504.71224793.33309231.4167142−3.193369913.475730323.2510209554.71225763.33311491.4167237−3.198085613.488562423.3032108Table 5
Comparison of invariants for the interaction of two solitary waves with results from [9] with h=0.1, Δt=0.025, A1=1, A2=0.5, (0≤x≤80).
tC1C2C3C1[9]C2[9]C3[9]04.71237323.33332531.41666434.71238843.33528901.416669754.71238613.33334821.41668524.71237183.33526351.4166486104.71239593.33336211.41669824.71238533.33528361.4166647154.71240653.33337851.41671414.71237563.33528941.4166772204.71242493.33341641.41675214.71237483.33530411.4166926254.71248993.33358321.41692384.71241733.33542781.4168363304.71276433.33335571.41776174.71264103.33594641.4176398354.71304743.33525001.41888494.71283533.33642471.4186746404.71248813.33363161.41716904.71239463.33559511.4170695454.71230023.33318781.41675804.71222733.33523641.4166637504.71224793.33309231.41671424.71215673.33511751.4165797554.71225763.33311491.41672374.71214003.33508471.4165527Interaction of two solitary waves at different times.
(a)(b)(c)(d)(e)(f)Interaction of two solitary waves at different times.
(a)(b)(c)(d)(e)(f)


## 4. Conclusion

In this paper, the cubic B-spline lumped Galerkin method has been successfully applied to obtain the numerical solution of the modified equal width wave equation. The efficiency of the method was tested on two test problems of wave propagation: the motion of a single solitary wave and the interaction of two solitary waves, and its accuracy was shown by calculating error normsL2 and L∞. It is clear that the error norms are adequately small and the invariants are satisfactorily constant in all computer run. The method can be also efficiently applied for solving a large number of physically important non-linear problems.



---

*Source: 587208-2012-01-16.xml*