On the Numerical Solution of One-Dimensional Nonlinear Nonhomogeneous Burgers ’ Equation

The nonlinear Burgers’ equation is a simple form of Navier-Stocks equation. The nonlinear nature of Burgers’ equation has been exploited as a useful prototype differential equation for modeling many phenomena. This paper proposes two meshfree methods for solving the one-dimensional nonlinear nonhomogeneous Burgers’ equation. These methods are based on the multiquadric (MQ) quasi-interpolation operator LW2 and direct and indirect radial basis function networks (RBFNs) schemes. In the present schemes, theTaylors series expansion is used to discretize the temporal derivative and the quasi-interpolation is used to approximate the solution function and its spatial derivatives. In order to show the efficiency of the present methods, several experiments are considered. Our numerical solutions are compared with the analytical solutions as well as the results of other numerical schemes. Furthermore, the stability analysis of the methods is surveyed. It can be easily seen that the proposed methods are efficient, robust, and reliable for solving Burgers’ equation.


Introduction
In this paper, we consider the one-dimensional nonlinear nonhomogeneous Burgers' equation: +   =   +  (, ) ,  ∈ Ω = [, ] ,  ⩾  0 , (1) with the initial condition, and the boundary conditions,  (, ) =  1 () , (, ) =  2 () , where (),  1 (), and  2 () are known functions,  is the positive parameter that related to the Reynolds number  = 1/, and (, ) is a known nonhomogeneous term.This equation was first derived from the hydrodynamics equations and used in surveying the laser generation of sound [1].Later on, it was applied to other physical phenomena such as wind forcing the buildup of water waves, electrohydrodynamic field in plasma physics, and design of feedback control [2][3][4].
When (, ) = 0, (1) is the well-known Burgers' equation: Burgers' equation in homogeneous form was first introduced by Bateman [5] who found its steady solutions, descriptive of certain viscous flows.It was later presented by Burgers as one of class of equations describing mathematical models of turbulence [6].In the context of gas dynamics, it was surveyed by Hopf [7] and Cole [8].The homogeneous Burgers' equation appears in various areas of applied mathematics and physics such as the phenomena of turbulence and supersonic flow, heat conduction, elasticity, and fusion [6][7][8][9][10].
From an analytical point of view, the nonhomogeneous form is poorly studied, the complete analytical solution being closely dependent on the form of the nonhomogeneous term.For example, Karabutov et al. [1] obtained the analytical solution of the nonhomogeneous Burgers' equation with (, ) =  sin(),  > 0, Ding et al. [11] studied the solution of (1) for the time-independent nonhomogeneous term (, ) = −, and Rao and Yadav [12] represented the solutions of the nonhomogeneous Burgers' equation for the nonhomogeneous term (, ) = /(2 + 1) 2 that  > 0 and  > 0 are constant.Recently, Moreau and vallée have obtained the analytical solution of the nonhomogeneous Burgers' equation with an elastic forcing term (, ) = −+ (),  ∈ R [3].
Up to now, various numerical methods are presented for the homogeneous Burgers' equation such as finite difference, finite element, boundary element, and collocation methods.For a survey of these methods refer to [13][14][15][16] and references cited therein.Among the methods that are mentioned above, the spatial domain where the partial differential governing equations are defined is often discretized into meshes.In these methods, the creation of suitable meshes is very essential for getting accurate results.However, the procedure of mesh generation consumes a lot of time and labor for some problems, especially for discontinuous and high gradient problems, for which these methods will become complicate.The root of these difficulties is the use of mesh in the formulation step.To avoid the mesh generation, recently, a kind of so-called meshfree or meshless method has extended quickly.In these methods, the scattered nodes are only used instead of meshing the domain of the problem.
For the last 20 years, the radial basis functions (RBFs) method was known as a powerful tool for the scattered data interpolation problem.The use of RBFs as a meshless process for the numerical solution of partial differential equations (PDEs) is based on the collocation scheme.The meshless methods based on RBFs were studied for approximating the solution of PDEs since initial development of Kansa's work (1990) [17].Kansa's method was extended to solve various ordinary and partial differential equations [18][19][20].In these works, the solution function is decomposed into RBFs and its derivatives are then arrived through differentiation that caused the reduction in convergence rate.In order to avoid this problem, Mai-Duy and Tran-Cong proposed an integrated MQ-RBFNs scheme for the approximation of function and its derivatives [21].Numerical experiments and theoretical analysis demonstrate that for solving PDEs integrated RBF (IRBF) procedure is more accurate in comparison with direct RBF (DRBF) procedure.Also, IRBF scheme is more stable than DRBF for a range of PDEs [21,22].
In both DRBF and IRBF schemes, one must resolve a linear system of equations at each time step.In the past decade, the other meshless method was introduced by using a MQ quasi-interpolation without solving a linear system of equations.MQ quasi-interpolation is constructed directly from linear combination of MQ-RBF and the approximated function.In 1992, Beaston and Powell [23] presented three univariate MQ quasi-interpolations named as L A , L B , and L C .Wu and Schaback [24] proposed the MQ quasiinterpolation L D and indicated that the scheme is shape preserving and convergent.Recently, Jiang et al. [25] have introduced a new multilevel univariate MQ quasi-interpolation approach with high approximation order compared with initial MQ quasi-interpolation scheme named as L W and L W 2 .This approach is based on inverse multiquadric (IMQ) RBF interpolation and Wu and Schaback's MQ quasi-interpolation operator L D .Chen and Wu applied initial MQ quasiinterpolation scheme for solving one-dimensional nonlinear homogeneous Burgers' equation [26].
In numerical solution of time dependent PDEs, such as Burgers' and Sine-Gordon equations, by using MQ quasiinterpolation scheme, there is a limitation for discretization of the temporal derivative.One has to use low order finite difference approximation for discretization of time derivatives because one does not solve any system of equations at each time step; otherwise one must solve a system of equations [27].Also, large number of nodes must be used for getting appropriate accuracy; see [26,28,29].
In this paper, we present two numerical methods by using MQ quasi-interpolation for the numerical solution of the nonhomogeneous Burgers' equation.In both of them, we use a two-order approximation for discretization of the time derivative.The main idea behind the discretization is to use more time derivatives in Taylor series expansion.This approach was demonstrated by Lax and Wendroff in finite difference [30] and used by Dag et al. for the homogeneous Burgers' equation [31].By using this discretization, we have to solve a linear system of equations at each time step that the size of the system is equivalent to the number of the centers in the spatial domain.Also, because the IRBFN scheme requires fewer centers in comparison with the DRBFN scheme, we apply MQ quasi-interpolation scheme in the indirect form in order not to encounter with large scale matrix.
The Jiang et al.MQ quasi-interpolation operator L W 2 is summation of two series that the second series coefficients are combined with first series coefficients.By giving relation between two series coefficients based on function values, we can convert it to a compact form based on one series and use it in direct and indirect forms for the numerical solution of PDEs.
The rest of present paper is organized as follows.A brief explanation of the MQ quasi-interpolation scheme is given in Section 2. Our numerical methods are applied on the nonlinear Burgers' equation in Section 3. In Section 4, the stability analysis of the methods is discussed.The results of several numerical experiments are reported in Section 5. Finally, some conclusions based on obtained results are drawn in Section 6.

The MQ Quasi-Interpolation Scheme
In this section, three univariate MQ quasi-interpolation schemes named as L D , L W , and L W 2 are described.Then, we describe our approach which converts operator L W 2 to the compact form.More details can be seen in [23][24][25].
For a given region Ω = [, ] and a finite set of distinct points, if we are supplied with a function  : [, ] → R, quasiinterpolation of  takes the form: where each function   () is a linear combination of the Hardy MQs basis function [32], and low order polynomials and  ∈ R + is a shape parameter.This formula requires the derivative values of  at the end points that is not convenient for practical purposes [23].Wu and Schaback [24] presented the univariate MQ quasiinterpolation operator L D that is defined as where Suppose that {   }  =1 is a smaller set from the given points {  }  =0 , where  is a positive integer satisfying  <  and 0 =  0 <  1 < ⋅ ⋅ ⋅ <  +1 = .Using the IMQ-RBF, the second derivative of () can be approximated by RBF interpolant    as where and  ∈ R + is a shape parameter.The coefficients {  }  =1 are uniquely determined by the interpolation condition Since ( 12) is solvable [33], so where By using the  and the coefficient  defined in (13), a function () is constructed in the form Now, the MQ quasi-interpolation operator L W by using L D defined by ( 8) and ( 9) on the data {(  , (  ))}  =1 with the shape parameter  is given by The shape parameters  and  should not be the same constant in (16).In (12), the value of      can be replaced by when the data's {(   , (   ))}  =1 are given.So, if    in ( 13) is replaced by then the quasi-interpolation operator defined by ( 15) and ( 16) is denoted by L W 2 The linear reproducing property and the high convergence rate of L W 2 were also studied in [25].
Therewith, ( 23) can be rewritten as Hence, the basis function ψ () are arrived as follows: By writing operator L W 2 in the compact form (23), we can use it in two indirect and direct forms for the numerical solution of PDEs.

The Numerical Methods
In this section, the numerical schemes are presented for solving the nonlinear Burgers' equation ( 1) by using the MQ quasi-interpolation L W 2 In our approach, the MQ quasiinterpolation approximates the solution function and the spatial derivatives of the differential equation and Taylor's series expansion is employed to approximate the temporal derivative similar to the work that Dag et al. did in [31].The MQ quasi-interpolation method is applied in direct and indirect forms.
We discretize the problem using the following Taylor's series expansion with step size Δ: where    =   (,   ) and   =  0 + Δ.Differentiating (1) with respect to time,    can be written as follows: where    =   (,   ).Using forward difference formula for the time derivative    in ( 27),    can be rewritten as Substituting ( 28) into (26) and using the expression achieved in (1), the following time discretized form of nonlinear Burgers' equation is yielded: 3.1.The Direct MQ Quasi-Interpolation Scheme.In this scheme, the unknown function   is approximated by using MQ quasi-interpolation scheme, and its spatial derivatives    and    are calculated by differentiating such closed form of quasi approximation as follows: where  ψ / = ψ  and  2 ψ / 2 =   .Now, replacing ( 30)-( 32) into ( 29) and applying collocation method yield where    = (  ,   ),    (  ) =   (  ,   ), ψ = ψ (  ), ψ  = ψ  (  ), and   =   (  ), whereas, according to (3), we have Substituting (34) and ( 35) into (33), wherein (33) generates a system of −1 linear equations in −1 unknown parameters  +1  .Equation ( 33) can be written in the matrix form where symbol * stands for component by component multiplication, Subsequently, (36) can be written as where , In order to make reduction in error, the obtained   from (38) is substituted in the right hand side of (30) that can be written as follows: and the obtained value is considered as   .Therefore, from (38) and (40), it yields that Hence, the unknown parameters   are specified from (41) instead of (38).

3.2.
The Indirect MQ Quasi-Interpolation Scheme.In indirect scheme, the highest order derivatives (second order in this paper) of the solution function are first approximated by (19), and their lower order derivatives and the solution function are then obtained by symbolic integration.Therefore,    can be approximated by MQ quasi-interpolation L W 2 on data {  } −1 =1 as follows: Now, integrating (42) yields Equations ( 42)-( 44) can be rewritten in the compact form as follows: where Similar to direct scheme, replacing (45) into ( 29) and ( 3) and applying collocation method lead to where Equations (47) generate a system of +1 linear equations in  + 1 unknown parameters  +1  .Similar to the direct quasi-interpolation scheme, (47) can be written in matrix form where, in this case, for  = 1, . . .,  − 1;  = 0, 1, . . .,  and for  = 0, ;  = 0, 1, . . .,  and Subsequently, (48) can be written as where and A 2 = A  + A  .From (42), it yields that Hence, the combination of ( 52) and ( 54) is given as

The Stability Analysis
In this section, the stability analysis from direct and indirect quasi-interpolation schemes is presented by using spectral radius of the amplification matrix similar to the work that Siraj-ul-Islam et al. did in [34].Let u be the exact and u * the numerical solution of (1); then the error vector  +1 = u +1 − u * +1 in the direct and indirect quasi-interpolation schemes can be written as where For the stability of the numerical schemes, we must have   → 0 as  → ∞; that is, (E 1 ) ⩽ 1, (E 2 ) ⩽ 1, which is the necessary and sufficient condition for the numerical schemes to be stable, where (E 1 ) and (E 2 ) denote the spectral radius of the amplification matrices E 1 and E 2 , respectively.Equations (56) can be written as Equations ( 57) can be written into the following forms by using the values of M 1 , N 1 , M 2 , and N 2 defined in (39) and (53): where The condition of stability will be satisfied if maximum eigenvalue of the matrix , and K depend on the shape parameter and the number of collocation points.Hence, the condition number and the spectral radius of the matrices E 1 and E 2 are dependent on the shape parameter and the number of collocation points.Since it is not possible to find explicit relationship among the spectral radius of the matrices and the shape parameter, this dependency is approximated numerically by keeping the number of collocation points fixed.

The Numerical Experiments
Five test experiments are studied to investigate the robustness and the accuracy of the proposed methods.The solution function of Burgers' equation is approximated by direct MQ quasi-interpolation (DMQQI) and indirect MQ quasiinterpolation (IMQQI) schemes and the results are compared with analytical solutions and the results in [18,26,31,35].The  ∞ and  2 error norms which are defined by are used to measure the accuracy.Also, the stability analysis of the methods is considered for first experiment.In all experiments, the shape parameter  is considered twice the shape parameter  and  is chosen twice .Also, the centers and the collocation points have been chosen as the same and equidistant.
The computations associated with the experiments discussed above were performed in Maple 16 on a PC with a CPU of 2.4 GHZ.
Experiment 1.In this experiment, we consider nonlinear Burgers' equation (1) with (, ) = 0 and the initial and the boundary conditions: The exact series solution of this experiment was given by Cole [8]: where Numerical results are presented for  = 0.1 and  = 0.01 with Δ = 0.001 and compared with the exact solutions and the results of the MQ quasi-interpolation scheme (MQQI; see [26]) and adaptive MQ scheme (AMQ; see [18]) for the cases of  = 0.1 and  = 0.01 in Tables 1 and 2, respectively.Also, the numerical solutions are compared with the results obtained by MQQI scheme, AMQ scheme, and Galerkin scheme [35] for  = 0.0001 in Table 3.In the cases  = 0.1 and  = 0.01, the shape parameter  is 0.815ℎ.In the case  = 0.0001, the parameter  is 2.78 × 10 −1 and 1.389 × 10 −4 for  = 36 and  = 72, respectively.The space-time graph of the estimated solution for  = 0.1 and  = 0.01 is presented in Figures 1 and  2.
Numerical comparison in these cases shows that the obtained results, particularly in IMQQI scheme, are in good agreement with the exact solutions and the results of the other schemes.
Relation between the spectral radius of the matrices E 1 and E 2 and the different values of the shape parameter  is shown in Table 4 by keeping the number of collocation points fixed.It is clear from Table 4 that if the values of shape parameter  are greater than the critical value  = 0.1 ( = 0.01), then the solution obtained from the IMQQI (DMQQI) method breaks down and hence the IMQQI and DMQQI methods become unstable.Therefore, the interval stability of IMQQI and DMQQI schemes is (0, 0.1) and (0.004, 0.01), respectively.
It can be seen from Table 4 that the schemes are very sensitive to the values of the shape parameter  and the interval stability of methods is a small interval.
Experiment 2. In this experiment, we consider the shock propagation solution of the homogeneous Burgers' equation [31] as a numerical experiment.This solution is given by The initial condition of the problem is obtained from (65) at time  = 1 and the boundary conditions in (3) can be obtained from the exact solution.Propagation of the shock is studied Table 1: Comparison of results with the exact solution and the results in [26] of  = 0.1 with Δ = 0.001 at  = 1 for different values of  of Experiment 1.    with Δ = 0.01 for  = 0.005 and  = 0.001.The shape parameter  is denoted by 6.0×10 −2 , 2.4×10 −2 , and 1.2×10 −2 for  = 20,  = 50, and  = 100, respectively.The  2 and  ∞ error norms are calculated in 240 and 2400 points for  = 0.005 and  = 0.001, respectively, and compared with the results of [31] in Table 5 at different times.The space-time graph of the estimated solution for  = 0.005 is showed in Figure 3.
Experiment 3. In this experiment, we study the fusion phenomenon of the two solitary waves of the homogeneous Burgers' equation.The fusion phenomenon happens when two or more solitons will fusion to one soliton at a specific time.In [9], Wang et al. studied the following Burgers' equation: They obtained the two-solitary-wave solution where  1 and  2 are constant.Let Hence, (66) converts to Burgers' equation form (1) wherein (, ) = 0 and  = 0.25.In this case, two-solitary-wave fusion happens at a specific time  = 0.Because we can show the fusion phenomenon, we consider an interval [−5, 5] for  and .For this purpose, we introduce a new time variable  =  + 5 and approximate the solution (, ) of (1) by using our schemes for  ∈ [0, 10].Then, we obtain (, ) for  ∈ [ −5, 5].The initial condition can be obtained from the exact solution at  = 0.The boundary conditions can be also taken from the exact solution.
Table 5: The comparison of  2 and  ∞ errors between the numerical results by using our schemes and the results of [31] with Δ = 0.01,  = 0.005, and  = 0.001 of Experiment 2. The  2 and  ∞ errors of IMQQI scheme are calculated in 1000 points for  1 = 1,  2 = −1, Δ = 0.001, and  = ℎ and listed in Table 6.From our numerical experiment whose results are not given here, we see that the accuracy of the DMQQI scheme is very bad in this experiment but we can clearly see from Table 6 that the results of the IMQQI scheme are in good agreement with the exact solutions.
The space-time graph of the estimated solution by using IMQQI is presented in Figure 4.
In this paper, we simulate this solution for  = 5 and  = 2 in  ∈ [−1, 1].The boundary conditions in (3) can be obtained from the exact solution (70).The  2 and  ∞ error norms are calculated in 1000 points for arbitrary , Δ = 0.01, and  = ℎ and listed in Table 7. Also, the space-time graph of the estimated solution by using IMQQI and DMQQI schemes is presented in Figure 5.
Table 7 shows that the accuracy of the DMQQI scheme is low even if the number of collocation points increases,     whereas the IMQQI scheme provides the good results with a small number of points.
The numerical results compared with the exact solutions for  = 1 with  = 20 and  = 0.1 with  = 100 in Tables 8 and 9, respectively.The numerical calculations are performed with Δ = 0.001 and  = ℎ.The space-time graph of the estimated solutions is also presented for  = 0.1 in Figure 6.

Conclusion
In this paper, two numerical schemes based on high accuracy MQ quasi-interpolation scheme and RBFs approximation schemes (IRBF and DRBF approximation schemes) have been presented for solving the nonlinear nonhomogeneous Burgers' equation.The accuracy of the methods can be increased by selecting the appropriate shape parameter.The choice of the shape parameter is still a pendent question.
The numerical results which are given in the previous section indicate that the performance of the methods specially IMQQI is in excellent agreement with the exact solutions.Tables 1-9 show that the IMQQI scheme is more accurate than DMQQI scheme as expected, and the interval stability of the IMQQI method is greater than the DMQQI method.Also, the IMQQI scheme required less nodes in comparison with the DMQQI scheme.We can even get good results with less number of points specially in IMQQI method but the results are bad at the ends of interval that we can improve it by using the knot method [37].Hence, we will not encounter with large scale matrix.Besides, we use equidistant points in our numerical experiments but our schemes can be used for the scattered points.
R 2 ] and maximum eigenvalue of the matrix E 2 = [I + (Δ/2)S 1 ] −1 [K + (Δ/2)S 2 ] are less than unity (in direct and indirect MQ quasi-  2 ,   1 ,   2 , and   denote the eigenvalues of the matrices R 1 , R 2 , S 1 , S 2 , and K, respectively.It is clear from (60) that the stability of the methods depends on the time step Δ and eigenvalues of the matrices   1 ,   2 ,   1 ,   2 , and   .The condition numbers and magnitude of the eigenvalues of the matrices R 1 , R 2 , S 1 , S 2

Table 2 :
[26]arison of results with the exact solution and the results in[26]of  = 0.01 with Δ = 0.001 at  = 1 for different values of  of Experiment 1.

Table 6 :
The  2 and  ∞ errors of the IMQQI scheme with Δ = 0.001 at different times of Experiment 3.

Table 7 :
The comparison of  2 and  ∞ errors between the numerical results of our schemes with Δ = 0.01 of Experiment 4.

Table 8 :
Comparison of numerical results with the exact solutions for  = 1 with  = 20 and Δ = 0.001 at  = 3 of Experiment 5.

Table 9 :
Comparison of numerical results with the exact solutions for  = 0.1 with  = 100 and Δ = 0.001 at  = 3 of Experiment 5.