    # A semi-discrete numerical method for convolution-type unidirectional wave equations

Numerical approximation of a general class of nonlinear unidirectional wave equations with a convolution-type nonlocality in space is considered. A semi-discrete numerical method based on both a uniform space discretization and the discrete convolution operator is introduced to solve the Cauchy problem. The method is proved to be uniformly convergent as the mesh size goes to zero. The order of convergence for the discretization error is linear or quadratic depending on the smoothness of the convolution kernel. The discrete problem defined on the whole spatial domain is then truncated to a finite domain. Restricting the problem to a finite domain introduces a localization error and it is proved that this localization error stays below a given threshold if the finite domain is large enough. For two particular kernel functions, the numerical examples concerning solitary wave solutions illustrate the expected accuracy of the method. Our class of nonlocal wave equations includes the Benjamin-Bona-Mahony equation as a special case and the present work is inspired by the previous work of Bona, Pritchard and Scott on numerical solution of the Benjamin-Bona-Mahony equation.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

In this paper, we propose a semi-discrete numerical approach based on a uniform spatial discretization and truncated discrete convolution sums for the computation of solutions to the Cauchy problem associated to the one-dimensional nonlocal nonlinear wave equation, which is a regularized conservation law,

 ut+(β∗f(u))x=0, (1.1)

with a general kernel function and the convolution integral

 (β∗v)(x)=∫Rβ(x−y)v(y)dy.

We prove error estimates showing the first-order or second-order convergence, in terms of the mesh size, depending on the smoothness of the kernel function. Also, our numerical experiments confirm the theoretical results.

Members of the class (1.1) arise as model equations in different contexts of physics, from shallow water waves to elastic deformation waves in dense lattices. For instance, in the case of the exponential kernel (which is the Green’s function for the differential operator where represents the partial derivative with respect to ) and with , (1.1) reduces to a generalized form of the Benjamin-Bona-Mahony (BBM) equation ,

 ut+ux−uxtt+(up+1)x=0, (1.2)

that has been widely used to model unidirectional surface water waves with small amplitudes and long wavelength. On the other hand, if the kernel function is chosen as the Green’s function of the differential operator ,

 β(x)=12√2e−|x|√2(cos(|x|√2)+sin(|x|√2)), (1.3)

and if , (1.1) reduces to a generalized form of the Rosenau equation 

 ut+ux+uxxxxt+(g(u))x=0 (1.4)

that has received much attention as a propagation model for weakly nonlinear long waves on one-dimensional dense crystal lattices within a quasi-continuum framework. It is worth to mention here that, in general, does not have to be the Green’s function of a differential operator. In other words, (1.1

) cannot always be transformed into a partial differential equation and those members might be called ”genuinely nonlocal”. Naturally, in such a case, standard finite-difference schemes will not be applicable to those equations. In this work we consider a numerical scheme based on truncated discrete convolution sums, which also solves such genuinely nonlocal equations.

We remark that the Fourier transform of the kernel function gives the exact dispersion relation between the phase velocity and wavenumber of infinitesimal waves in linearized theory. In other words, general dispersive properties of waves are represented by the kernel function. Obviously the waves are nondispersive when the kernel function is the Dirac delta function. That is, (

1.1) can be viewed as a regularization of the hyperbolic conservation law , where the convolution integral plugged into the conservation law is the only source of dispersive effects. Naturally, our motivation in developing the present numerical scheme also stems from the need to understand the interaction between nonlinearity and nonlocal dispersion.

Recently, two different approaches have been proposed in [3, 4] to solve numerically the nonlocal nonlinear bidirectional wave equation . In  the authors have developed a semi-discrete pseudospectral Fourier method and they have proved the convergence of the method for a general kernel function. By pointing out that, in most cases, the kernel function is given in physical space rather than Fourier space, in  the present authors have developed a semi-discrete scheme based on spatial discretization, that can be directly applied to the bidirectional wave equation with a general arbitrary kernel function. They have proved a semi-discrete error estimate for the scheme and have investigated, through numerical experiments, the relationship between the blow-up time and the kernel function for the solutions blowing-up in finite time. These motivate us to apply a similar approach to the initial-value problem of (1.1) and to develop a convergent semi-discrete scheme that can be directly applied to (1.1). To the best of our knowledge, no efforts have been made yet to solve numerically (1.1) with a general kernel function.

As in  our strategy in obtaining the discrete problem is to transfer the spatial derivative in (1.1) to the kernel function and to discretize the convolution integral on a uniform grid. Thus, the spatial discrete derivatives of

do not appear in the resulting discrete problem. Depending on the smoothness of the kernel function, the two error estimates, corresponding to the first- and second-order accuracy in terms of the mesh size, are established for the spatially discretized solution. If the exact solution decreases fast enough in space, a truncated discrete model with a finite number of degrees of freedom (that is, a finite number of grid points) can be used and, in such a case, the number of modes depends on the accuracy desired. Of course, this is another source of error in the numerical simulations and it depends on the decay behavior of the exact solution. Following the idea in

 we are able to prove a decay estimate for the exact solution under certain conditions on the kernel function. We also address the above issues in two model problems; propagation of solitary waves for both the BBM equation and the Rosenau equation.

A numerical scheme based on the discretization of an integral representation of the solution was used in  to solve the BBM equation which is a member of the class (1.1). The starting point of our numerical method is similar to that in . As it was already observed in , an advantage of this direct approach is that a further time discretization will not involve any stability issues regarding spatial mesh size. This is due to the fact that our approach does not involve any spatial derivatives of the unknown function.

The paper is structured as follows. In Section 2 we focus on the continuous Cauchy problem for (1.1). In Section 3, the semi-discrete problem obtained by discretizing in space is presented and a short proof of the local well-posedness theorem is given. In Section 4 we investigate the convergence of the discretization error with respect to mesh size; we prove that the convergence rate is linear or quadratic depending on the smoothness of . Section 5 is devoted to analyzing the key properties of the truncation error arising when we consider only a finite number of grid points. In Section 6 we carry out a set of numerical experiments for two specific kernels to illustrate the theoretical results.

The notation used in the present paper is as follows. is the () norm of on , is the -based Sobolev space with the norm and is the usual -based Sobolev space of index on . denotes a generic positive constant. For a real number , the symbol denotes the largest integer less than or equal to .

## 2 The Continuous Cauchy Problem

We consider the Cauchy problem

 ut+(β∗f(u))x=0, \ \ \ x∈R,% \ \ t>0, (2.1) u(x,0)=φ(x), \ \ \         x∈R. (2.2)

We assume that is sufficiently smooth with and that the kernel satisfies:

1. ,

2. is a finite Borel measure on .

We note that Condition 2 above also includes the more regular case in which case . The following theorem deals with the local well-posedness of (2.1)-(2.2).

###### Theorem 2.1

Suppose that satisfies Assumptions 1 and 2. Let , with . For a given , there is some so that the initial-value problem (2.1)-(2.2) is locally well-posed with solution .

The proof of Theorem 2.1 follows from Picard’s theorem for Banach space-valued ODEs. The nonlinear term is locally Lipschitz on . Moreover, the conditions on imply that the term maps onto itself. Hence, (2.1) is an -valued ODE.

We will later use the following estimate on the nonlinear term .

###### Lemma 2.2

Let with . Then for any , we have . Moreover there is some constant depending on such that for all with

 ∥f(u)∥Hs≤C(M)∥u∥Hs .

We emphasize that the bound in the lemma depends only on the -norm of . This in turn allows us to control the -norm of the solution by its -norm. In particular, finite time blow-up of solutions is independent of regularity and is controlled only by the -norm of .

###### Lemma 2.3

Suppose the conditions of Theorem 2.1 are satisfied and is the solution of (2.1)-(2.2). Then

 ∥u(t)∥Hs≤∥φ∥HseCt,    0≤t≤T, (2.3)

where depends on and .

Proof. Integrating (2.1) with respect to time, we get

 u(x,t)=φ(x)−∫t0(β′∗f(u))(x,τ)dτ. (2.4)

Using Young’s inequality and Lemma 2.2, we obtain

 ∥u(t)∥Hs≤∥φ∥Hs+C(M,β)∫t0∥u(τ)∥Hsdτ (2.5)

for all . By Gronwall’s lemma this gives (2.3).

###### Remark 2.4

It is well-known that, under suitable convexity assumptions, the hyperbolic conservation law leads to shock formation in a finite time even for smooth initial conditions. On the other hand, according to Lemma 2.3, for any solution of (1.1) the derivatives will stay bounded as long as stays bounded in the -norm. In other words, the regularization prevents a shock formation.

## 3 The Discrete Problem

In this section we first give two lemmas for error estimates of discretizations of integrals and derivatives on an infinite uniform grid, respectively, and then introduce the discrete problem associated to (2.1)-(2.2).

### 3.1 Discretization and Preliminary Lemmas

Consider doubly infinite sequences of real numbers with (where denotes the set of integers). For a fixed and , the space is defined as

 lph(Z)={(wi):wi∈R,  ∥w∥plph=∞∑i=−∞h|wi|p}.

The space with the sup-norm is a Banach space. The discrete convolution operation denoted by the symbol transforms two sequences and into a new sequence:

 (w∗v)i=∑jhwi−jvj (3.1)

(henceforth, we use to denote summation over all ). Young’s inequality for convolution integrals state that for , and for .

Consider a function of one variable defined on . We then introduce a uniform partition of the real line with the mesh size and with the grid points , . Let the restriction operator be . We will henceforth use the abbreviations and for and , respectively.

The following lemma gives the error bounds for the discrete approximations of the integral over depending on the smoothness of the integrand. We note that the two cases correspond to the rectangular and trapezoidal approximations, respectively.

###### Lemma 3.1
1. Let and be a finite measure on . Then

 ∣∣ ∣∣∫Rw(x)dx−∑ihw(xi)∣∣ ∣∣≤h|μ|(R). (3.2)
2. Let and be a finite measure on . Then

 ∣∣ ∣∣∫Rw(x)dx−∑ihw(xi)∣∣ ∣∣≤h2|ν|(R). (3.3)

The following lemma handles the estimates for the discrete approximation of the first derivative , whose proof follows more or less standard lines.

###### Lemma 3.2

Let and . Let be the discrete derivative operator defined by the central differences

 (Dw)i=12h(wi+1−wi−1),    i∈Z. (3.4)
1. If , then

 ∥Dw−w′∥l∞≤h2∥w′′∥L∞. (3.5)
2. If , then

 ∥Dw−w′∥l∞≤h26∥w′′′∥L∞. (3.6)

We refer the reader to  for the proofs of the two lemmas above.

### 3.2 The Semi-Discrete Problem

In order to get the semi-discrete problem associated with (2.1)-(2.2) we discretize them in space with a fixed mesh size . Thus, the discretized form of the nonlocal wave equation (2.1) becomes

 dvdt=−D(βh∗f(v)) (3.7)

with the notation and . The identity allows us to transfer the discrete derivative on the kernel. So in order to prove the local well-posedness theorem for the semi-discrete problem, we need merely prove the following lemma which estimates the discrete derivative of the restriction of the kernel.

###### Remark 3.3

We note that (3.7) involves point values of . When satisfies Assumptions 1 and 2 one should pay attention to how the point values in (3.7) are defined. To clarify this issue we will assume throughout that

 β(x)=∫(−∞,x]dμ. (3.8)
###### Lemma 3.4

Let be a finite measure on . Then and

Proof. The assumption (3.8) allows one to write

 h(Dβh)i =12(β(xi+1)−β(xi−1)) (3.9) =12∫(xi−1,xi+1]dμ≤12|μ|((xi−1,xi+1]), (3.10)

from which we deduce the estimate

Let be a locally Lipschitz function with . The map is locally Lipschitz on . Moreover, by Lemma 3.4, the map

 v⟶D(βh∗f(v))

is also locally Lipschitz on . By Picard’s theorem on Banach spaces, this implies the local well-posedness of the initial-value problem for (3.7).

###### Theorem 3.5

Let be a locally Lipschitz function with . Then the initial-value problem for (3.7) is locally well-posed for initial data in . Moreover there exists some maximal time so that the problem has unique solution . The maximal time , if finite, is determined by the blow-up condition

 limsupt→T−h∥v(t)∥l∞=∞. (3.11)

## 4 Discretization Error

Suppose that the function with sufficiently large is the unique solution of the continuous problem (2.1)-(2.2). The discretizations of the initial data on the uniform infinite grid will be denoted by . Let be the unique solution of the discrete problem based on (3.7) and the initial data . The aim of this section is to estimate the discretization error defined as . Depending on the conditions imposed on the kernel function , we provide the two different theorems establishing the first- and second-order convergence in . The proofs follow similar lines as the corresponding ones in .

###### Theorem 4.1

Suppose that satisfies Assumptions 1 and 2. Let , with . Let be the solution of the initial-value problem (2.1)-(2.2) with . Similarly, let be the solution of (3.7) with initial data . Let . Then there is some so that for , the maximal existence time of is at least and

 ∥u(t)−uh(t)∥l∞=O(h) (4.1)

for all .

Proof. We first let . Since , by continuity there is some maximal time such that for all . Moreover, by the maximality condition either or .  At the point , (2.1) becomes

 ut(xi,t)+(β∗f(u))x(xi,t)=0.

Recalling that , this becomes . A residual term arises from the discretization of (2.1):

 (4.2)

where

The th entry of satisfies

 (Fh)i =(F1h)i+(F2h)i, (4.3)

where the variable is suppressed for brevity. We start with the term . Replacing by for convenience, we have

 (F1h)i=(βh∗g′)i−(β∗g′)(xi)=∑jhβ(xi−xj)g′(xj)−∫Rβ(xi−y)g′(y)dy. (4.4)

Since and be a finite measure on , by (3.2) of Lemma 3.1 we have

 ∣∣(F1h)i∣∣≤h|˜μ|(R)

where with . When , we have

 ∥r′∥L1≤∥β′∥L1∥g′∥L∞+∥β∥L1∥g′′∥L∞.

Since and , are bounded for , we have

 |˜μ|(R) ≤|μ|(R)∥g′∥L∞+∥β∥L1∥g′′∥L∞, (4.5) ≤C(|μ|(R)+∥β∥L1)∥u∥Hs (4.6)

where we have used Lemma 2.2. Thus

 |(F1h)i|≤Ch(|μ|(R)+∥β∥L1)∥u∥Hs≤Ch∥u∥Hs.

For the second term , again with and we have

 ∣∣(F2h)i∣∣ =∣∣(βh∗(Dg−Rg′))i∣∣≤C∥Dg−g′∥l∞, ≤Ch∥g′′∥L∞≤Ch∥g∥Hs≤Ch∥u∥Hs,

where Lemma 2.2 and (3.5) of Lemma 3.2 are used. Combining the estimates for and , we obtain

 ∥Fh(t)∥l∞≤Ch∥u(t)∥Hs.

We now let be the error term. Then, from (3.7) and (4.2) we have

 de(t)dt=−Dβh∗(f(u)−f(uh))+Fh,    e(0)=0.

This implies

 e(t)=∫t0(−Dβh∗(f(u)−f(uh))+Fh)dτ. (4.7)

By noting that

 ≤∥Dβh∥l1h∥f(u)−f(uh)∥l∞ ≤C∥u−uh∥l∞≤C∥e(t)∥l∞,

it follows from (4.7) that, for ,

 ∥e(t)∥l∞ ≤sup0≤t≤T∥Fh(t)∥l∞∫t0dτ+C∫t0∥e(τ)∥l∞dτ, ≤C(ht+∫t0∥e(τ)∥l∞dτ). (4.8)

Then, by Gronwall’s inequality,

 ∥e(t)∥l∞≤ChTeCT.

We observe that the constant depends on the bounds , and . We note that, by Lemma 2.3, depends on and . The above inequality, in particular, implies that for sufficiently small . Then we have showing that . From the above estimate we get (4.1).

###### Theorem 4.2

Let , with . Let and be a finite measure on . Let be the solution of the initial-value problem (2.1)-(2.2) with . Similarly, let be the solution of (3.7) with initial data . Let . Then there is some so that for , the maximal existence time of is at least and

 ∥u(t)−uh(t)∥l∞=O(h2) (4.9)

for all .

Proof. The proof follows that of Theorem 4.1 very closely. So here we use the same notation and provide only a brief outline. The main observation needed is that the only place where the proofs differ is the estimate for in (4.4). Since and be a finite measure on , by the second part of Lemma 3.1 we have

 ∣∣(F1h)i∣∣≤h2|˜ν|(R)

where with . Formally we can write

 r′′(y)=β′′(xi−y)g′(y)−2β′(xi−y)g′′(y)+β(xi−y)g′′′(y).

Noting that and , , are bounded for and using Lemma 2.2, we get

 |˜ν|(R) ≤C(|ν|(R)∥g′∥L∞+2∥β′∥L1∥g′′∥L∞+∥β∥L1∥g′′′∥L∞ ≤C(|ν|(R)+2∥β∥W1,1)∥u∥Hs (4.10)

so that

For the second term , again with and we have

 ∣∣(F2h)i∣∣ =∣∣(βh∗(Dg−Rg′))i∣∣≤C∥Dg−g′∥l∞ ≤Ch2∥g′′′∥L∞≤Ch2∥g∥Hs≤Ch2∥u∥Hs,

where (3.6) of Lemma 3.2 is used. Using the estimates for and in (4.3), we obtain

 ∥Fh(t)∥l∞≤Ch2∥u∥Hs.

Now (4.8) takes the following form

 ∥e(t)∥l∞≤C(h2t+∫t0∥e(τ)∥l∞dτ) (4.11)

for . Then, by Gronwall’s inequality,

 ∥e(t)∥l∞≤Ch2TeCT.

Now the constant depends on the bounds , and . The rest of the proof follows the same lines as the proof of Theorem 4.1.

## 5 The Truncated Problem and a Decay Estimate

### 5.1 The Truncated Problem

In practical computations, one needs to truncate both the infinite series in (3.1) at a finite and the infinite system of equations in (3.7) to the system of equations. After truncating, we obtain from (3.7) the finite-dimensional system

 dvNidt=−N∑j=−NhDβ(xi−xj)f(vNj),% \ \ \ \ \ \ −N≤i≤N (5.1)

where

are the components of a vector valued function

with finite dimension . In this section we estimate the truncation error resulting from considering (5.1) instead of (3.7) and give a decay estimate for solutions to the initial-value problem (2.1)-(2.2) with certain kernel functions.

We start by rewriting (5.1) in vector form

 dvNdt=−BNf(vN),

where denotes a matrix with rows and columns, whose typical element is . Using the notation and assuming that is a finite measure on we get from Lemma 3.4

 ∥BNw∥l∞≤|μ|(R)∥w∥l∞.

As long as is a bounded and smooth function, there exists a unique solution of the initial-value problem defined for (5.1) over an interval . Moreover the blow-up condition

 limsupt→(TN)−∥vN(t)∥l∞=∞ (5.2)

is compatible with (3.11) in the discrete problem. Consider the projection of the solution of the semi-discrete initial-value problem associated with (3.7) onto , defined by

 TNv=(v−N,v−N+1,…,v0,…,vN−1,vN)

with the truncation operator . Our goal is to estimate the truncation error and show that, for sufficiently large , approximates the solution of the continuous problem (2.1)-(2.2).

###### Theorem 5.1

Let be the solution of (3.7) with initial value and let

 δ=sup{|vi(t)|:t∈[0,T],|i|>N}  and  ϵ(δ)=max|z|≤δ|f(z)|.

Then for sufficiently small , the solution  of (5.1) with initial value exists for times and

 ∣∣vNi(t)−vi(t)∣∣≤Cϵ(δ), \ t∈[0,T],

for all .

Proof. We follow the approach in the proof of Theorem 4.1. Taking the components with of (3.7) we have

 dvidt =−∞∑j=−∞hDβ(xi−xj)f(vj) =−N∑j=−NhDβ(xi−xj)f(vj)+FNi

with the residual term

 FNi=−∑|j|>NhDβ(xi−xj)f(vj).

Then satisfies the system

 dTNvdt=−BNf(TNv)+FN

with the residual term . Estimating the residual term we get

 ∣∣FNi∣∣≤|μ|(R)sup|j|>N∣∣f(vj)∣∣≤Cϵ(δ).

We set . Since , by continuity of the solution of the truncated problem there is some maximal time such that we have for all . By the maximality condition either or . We define the error term . Then

 d˜e(t)dt=−BN(f(TNv)−f(vN))+FN,    ˜e(0)=0,

so

 ˜e(t)=∫t0(−BN(f(TNv)−f(vN))+FN)dτ.

Then

 ∥˜e(t)∥l∞≤|μ|(R)∫t0∥∥(f(TNv)−f(vN))(τ)∥∥l∞dτ+∫t0∥FN(τ)∥l∞dτ.

But

 ∥∥f(TNv)−f(vN)∥∥l∞≤C∥∥TNv−vN∥∥l∞.

Putting together, we have

 ∥˜e(t)∥l∞≤CTϵ(δ)+C∫t0∥˜e(τ)∥l∞dτ,

and by Gronwall’s inequality,

 ∥˜e(t)∥l∞≤Cϵ(δ)TeCT.

This, in particular, implies that there is some such that for all we have . Then we have showing that and this completes the proof.

By combining Theorem 5.1 with Theorems 4.1 and 4.2, respectively, we now state our main results through the following two theorems.

###### Theorem 5.2

Let , with . Let and be a finite measure on . Let be the solution of the initial-value problem (2.1)-(2.2) with . Then for sufficiently small and , there is an so that the solution  of (5.1) with initial values , exists for times and

 ∣∣u(ih,t)−(uNh)i(t)∣∣=O(h+ϵ), \ t∈[0,T] (5.3)

for all .

###### Theorem 5.3

Let , with . Let and