Numerical Solutions of Coupled Burgers′ Equations

In this article, two new modified variational iteration algorithms are investigated for the numerical solution of coupled Burgers′ equations. These modifications are made with the help of auxiliary parameters to speed up the convergence rate of the series solutions. Three numerical test problems are given to judge the behavior of the modified algorithms, and error norms are used to evaluate the accuracy of the method. Numerical simulations are carried out for different values of parameters. The results are also compared with the existing methods in the literature.


Introduction
In recent years, coupled partial differential equations have been employed in various fields of engineering and applied sciences.Coupled Burgers' equations are coupled partial differential equations (PDEs), and describe the approximation theory of flow through a shock wave traveling in a viscous fluid.
In this paper, the following type of coupled Burger equations will be investigated: where α 1 and α 2 are positive and non-zero viscosity parameters, β 1,1 , β 1,2 , β 2,2 and β 2,1 are constants depend on the Stokes velocity.The initial and boundary conditions for coupled Burgers' Equation (1) are the following: Many authors have solved different types of coupled Burgers equations by different techniques such as the Fourier pseudospectral method [1] applied by Rashid and Ismail for one-dimensional coupled Burgers equations.A mesh-free interpolation method was employed by Islam et al. [2], where the radial basis function (RBF) collocation method has been formulated for the numerical treatment of coupled Burger equations and other nonlinear PDEs.Khater et al. [3] used the Chebyshev spectral collocation method, which has been developed using Chebyshev polynomials, and the Runge-Kutta 4th order method (RK4) method.Kya [4] used the decomposition method for the solution of a viscous coupled Burger equation, and obtained solutions in the form of a convergent power series.Mittal and Arora [5] proposed a scheme known as the Lai cubic B-spline collocation scheme for the solution of coupled viscous Burger equations, where the authors used a crank Nicholson scheme and cubic B-spline functions for time integration and space integration, respectively.Lai and Ma [6] proposed a new lattice Boltzmann model for the solution of coupled Burger equations, after selecting the distribution functions, Chapman-Enskog expansion was employed for the solution of coupled Burger equations.The Chebyshev-Legendre pseudospectral method [7] has been utilized by Rashid et al. for coupled viscous Burgers (VB) equations, where the leapfrog scheme and Chebyshev-Legendre Pseudo-Spectral method (CLPS) method were used for the time direction and space direction, respectively.Kumar and Pandit [8] used a composite numerical scheme for the numerical evaluation of coupled Burger equations, where the scheme was developed based on Haar wavelets and finite difference.Mohammadi and Mokhtari [9] used a reproducing Kernal method for an analytical solution in the form series of systems of Burger equations.At last, in 2019, Bak et al. [10] developed a new approach named a semi-Lagrangian approach for the numerical solution of coupled Burger equations.We compare our results with those of [10], and show the applicability of our proposed method to handle such problems as those that arise in applied science and engineering.
This paper aims to solve three types of coupled Burgers equations by employing variational iteration algorithm-II and one of the examples to be solved by modified variational iteration algorithm-I.The organization of the rest of the paper is as follows; in Section 2, we elucidate the variational iteration algorithm-II, in the next Section 3, the semi-numerical method is applied to three test problems, and a comparison is made with some other methods; lastly, some conclusions are drawn in the last Section 4.

Modified Variational Iteration Algorithm-II
Consider the nonlinear diffusion equation: where the terms L[w(x, t)] and N[w(x, t)] are linear and nonlinear terms in that order, while c(x, t) is the source term.For a given w 0 (x, t), the approximate solution w n+1 (x, t) of Equation ( 4) can be obtained as below: In Equation ( 5), λ is a parameter, which is known as the Lagrange multiplier [11].This Lagrange multiplier is obtained by variational theory, where w n is a term being restricted, which in turn gives δ w n = 0 and gives the following Lagrange multipliers: and in general, the Lagrange multiplier for m ≥ 1: Putting the value of the Lagrange multiplier from Equation (6) in the correctional function shown in Equation ( 5), we get the following iterative formula: This is known as variational iteration algorithm-I [12][13][14], which is an advance improvement of the common Lagrange multiplier method [11].These days, this technique [15][16][17][18][19][20][21][22][23] has been set up for offering a solution for a more extensive scope of problems, developing in several fields of pure and applied sciences.PDEs extensively arise in various physical applications such as propagation and the scattering of waves, magnetohydrodynamic flow through pipes, computational fluid dynamics, magnetic resonance imaging, the phenomena of turbulence and supersonic flow, the flow of a shock wave traveling in a viscous fluid, acoustic transmission, traffic and aerofoil flow theory, and the proposed technique has the ability to investigate these types problems effectively.
Inserting an auxiliary parameter in the variational iteration method and presenting a new modification of it, as below: we call it modified variational iteration algorithm-I.We will apply this method for a test problem for verifying its accuracy and compactness.According to He et al. [21], we can construct a more concise iteration formulation, which can be written as: Equation ( 9) is called variational iteration algorithm-II; we further modify it by coupling an auxiliary term, and named it modified variational iteration algorithm-II, which can be written as: The iterative algorithm does not contain unsure constants aside from an auxiliary parameter h, which is utilized to control the convergence of the obtained solution ideally by limiting the norm 2 of the residual error over the space of the given problem.The ideal decision of this h improves the precision and proficiency of the algorithm.In the wake of presenting h, Equation (10) takes the structure: w 0 (ζ) is a proper initial approximation, Now, we may start the procedures with the selective function w 0 (x, t) and use the above iteration structures to get the approximate solutions.The iteration structure shown in Equation (11) will give several approximations of w(x, t), and the exact solutions are obtained at the limit of the resulting successive approximations, i.e.: This algorithm is named MVIA-II.We utilize this modified algorithm for the solution of coupled nonlinear partial differential equations.

Numerical Assessments
In this section, we discuss three test problems of different types of coupled Burgers' equations to check the accuracy of the proposed method.We assess the accuracy of the method by taking different values of parameters, and the obtained results are very encouraging and significant, while the error norms are calculated and compared with the error of the other methods available in the literature.For numerical computation, Maple 16 and MATLAB R2015a are used on a Dell core i3, 1.90-GHz PC with 4 GB of RAM.

Test Problem 1
Consider the coupled Burger s Equation (1) with The equation becomes: The exact solution was given by [10]: where β 1,2 and β 2,1 are nonzero arbitrary constants, k = 0.05 and B = k . In this example, we investigate the proficiency and accuracy of the proposed algorithm by changing the values of arbitrary constants β 1,2 and β 2,1 , whose different magnitudes from each other produce the changed performance of u and v.The approximate numerical solutions are calculated with different parameters (β 1,2 , β 2,1 ) = (0.1, 0.3) by varying the t from t = 0.5 to 10.0.Constructing the correction function for the Equation ( 13) as: The values of λ 1 (ζ) and λ 2 (ζ) may be obtained most positively by the variational principle [24,25].We obtain the estimation of λ First, we use the modified variational iteration algorithm-I; then, we will use the modified algorithm-II for the solution of this test problem.Utilizing the estimation of λ 1 (ζ) and λ 2 (ζ) in Equations ( 15) and ( 16) results in the underneath iterative structure: Introducing with a proper initial guess: one can get the beneath different approximations by utilizing the recurrence relations shown in Equations ( 17) and ( 18) and β 1,2 = 0.1, β 2,1 = 0.3: . . .we stop the procedure at the third-order approximation.In order to find an appropriate value of h for the approximate solutions, a residual function can be defined as: The square of the r 3a (x, t, h) and r 3b (x, t, h) functions with respect to the auxiliary term h in the domain ( For finding an optimum value of the auxiliary term h, the lowest point of the error of norm 2 of Function ( 19) is selected.The lowest point of e 3a (h), as h = 0.994178085892076 and e 3b (h) as h = 0.996482304232335 is gained by means of Maple mathematical software.By replacing the value of the auxiliary term in u 3 (x, t, h), the absolute error of the fourth-order approximation of the proposed algorithm decreases remarkably.Now, we employ the modified variational iteration algorithm-II for the solution of test problem 1. Utilizing this estimation of λ(ζ) in Equation ( 13) results in the underneath iterative structure: Introducing with a proper initial guess: one can get the beneath different approximations by utilizing the recurrence relations shown in Equations ( 21) and ( 22) for β 1,2 = 0.1 and β 2,1 = 0.3: we stop the procedure at the third-order approximation.In order to find an appropriate value of h for the approximate solutions, a residual function can be defined as: The square of the r 3a (x, t, h) and r 3b (x, t, h) functions with respect to the auxiliary term h in the domain (x, t) [0, 1] × [0, 1] is: For finding an optimum value of the auxiliary term h, the lowest point of the error of norm 2 of Function ( 24) is selected.The lowest point of e 3a (h) as h = 1.00000000005237, and e 3b (h) as h = 0.999999981901321 , are gained by means of Maple mathematical software.By replacing the value of the auxiliary term in u 4 (x, t, h), the absolute error of the fourth-order approximation of the proposed algorithm decreases remarkably.A comparison is presented in the following tables with the Fourier pseudospectral method [1], meshfree interpolation method [2], Chebyshev spectral collocation method [3], decomposition method [4], cubic B-spline collocation scheme [5], new lattice Boltzmann model [6], and semi-Lagrangian approach [10].
In Tables 1 and 2, the L 2 and L ∞ error norms are given for our proposed algorithm MVIA-II, MVIA-I, and finite difference method (FDM), as well as the lattice Boltzmann model (LBM) for u and v, respectively, while Table 3 shows the comparison of errors of different methods for showing the accuracy of our proposed method.We can see that MVIA-II produces more accurate results than both FDM and LBM.The L 2 and L ∞ error norms for different values of t are reported in Tables 2 and 3, which show that the proposed method gives more accurate results than the previous methods that have been used before.Figures 1 and 2 display the behavior of exact and numerical solutions for u and v by MVIA-II. Figure 3 presents the absolute error graphs for different values of t for u and v, while Figure 4 is devoted to a plot solution at different time levels computed by MVIA-II.It can be seen that the method is uniformly convergent and gives accurate results.It is clear from the figures and tables that the proposed method can handle the coupled Burgers' equations very accurately.In Tables 1 and 2, the  and  error norms are given for our proposed algorithm MVIA-II, MVIA-I, and finite difference method (FDM), as well as the lattice Boltzmann model (LBM) for  and , respectively, while Table 3 shows the comparison of errors of different methods for showing the accuracy of our proposed method.We can see that MVIA-II produces more accurate results than both FDM and LBM.The  and  error norms for different values of t are reported in Tables 2  and 3, which show that the proposed method gives more accurate results than the previous methods that have been used before.Figures 1 and 2 display the behavior of exact and numerical solutions for  and  by MVIA-II. Figure 3 presents the absolute error graphs for different values of  for  and , while Figure 4 is devoted to a plot solution at different time levels computed by MVIA-II.It can be seen that the method is uniformly convergent and gives accurate results.It is clear from the figures and tables that the proposed method can handle the coupled Burgers' equations very accurately.

Test Problem 2
Consider the coupled Burgers' Equation ( 1) with The exact solution was given by [10]: where k is a non-zero arbitrary constant.We noticed that when the value of k increases, then the traveling waves become faster as time goes on, and the waveforms become increasingly steeper.Constructing the correction function for Equation (25) as: The values of λ 1 (ζ) and λ 2 (ζ) may be obtained most positively by the variational principle.We obtain the estimation of Utilizing this estimation of λ(ζ) in Equations ( 27) and (28) results in the underneath iterative structure: Introducing with a proper initial guess: one can get beneath different iterations by utilizing the recurrence relations shown in Equations ( 29) and (30), for k = 0.1: we stop the procedure at u 4 (x, t, h).In order to find an appropriate value of h for the approximate solutions, a residual function can be defined as: The square of the r 4 (x, t, h) function with respect to the auxiliary term h in the domain (x, t) [0, 1] × [0, 1] is: For finding an optimum value of the auxiliary term h, the lowest point of the error of norm 2 of Function (31) is selected.The lowest point of e 4 (h), as h = 1.00000000005237, is gained by means of Maple mathematical software.By replacing the value of the auxiliary term in u 4 (x, t, h), the absolute error of the fourth-order approximation of the proposed algorithm decreases remarkably.
In Tables 4 and 5, L 2 and L ∞ error norms are given for our proposed algorithm MVIA-II, and other methods available in the literature.We can see that MVIA-II produces more accurate results than FDM, LBM, and CBSLM.The L 2 and L ∞ error norms for different values of t are reported in Table 5, which shows that the proposed method gives more accurate results than previous methods used.Figure 5 displays the behavior of the exact and numerical solutions, while Figure 6 is devoted to a plot absolute error graph for different values of t.The method is uniformly convergent and gives accurate results.It is clear from the figures and tables that the proposed method can handle the coupled Burgers' equations very accurately.In Tables 4 and 5,  and  error norms are given for our proposed algorithm MVIA-II, and other methods available in the literature.We can see that MVIA-II produces more accurate results than FDM, LBM, and CBSLM.The  and  error norms for different values of  are reported in  In Tables 4 and 5,  and  error norms are given for our proposed algorithm MVIA-II, and other methods available in the literature.We can see that MVIA-II produces more accurate results than FDM, LBM, and CBSLM.The  and  error norms for different values of  are reported in Table 5, which shows that the proposed method gives more accurate results than previous methods used.Figure 5 displays the behavior of the exact and numerical solutions, while Figure 6 is devoted

Test Problem 3
Consider the coupled Burger Equation (1) with α The exact solution was given by [10]: where β 1,2 and β 2,1 are non-zero arbitrary constants.Constructing the correction function for Equation (33) as: The values of λ 1 (ζ) and λ 2 (ζ) may be obtained most positively by the variational principle.We obtain the estimation of λ 1 (ζ) and λ Utilizing this estimation of λ(ζ) in Equations ( 35) and (36) results in the underneath iterative structure: Introducing with a proper initial guess: one can get the beneath different approximations by utilizing the recurrence relations shown in Equations ( 37) and (38): The square of the r 11a (x, t, h) and r 11b (x, t, h) functions with respect to the auxiliary term h in the domain x ∈ [−π, π] is:  For finding an optimum value of the auxiliary term h, the lowest point of the error of norm 2 of Function (40) is selected.The lowest point of e 11 (h), as h = 1.00000000005237, is gained by means of Maple mathematical software.By replacing the value of the auxiliary term in u 11 (x, t, h), the absolute error of the 11th-order approximation of the proposed algorithm decreases remarkably.
In Table 6, L 2 and L ∞ error norms are given for our proposed algorithm MVIA-II and other numerical methods available in the literature.We can see that MVIA-II produces more accurate results than the others.Figure 7 displays the comparison of exact and numerical solutions of u and v for t = 1. Figure 8 is devoted to showing the behavior of exact and approximate solutions for u and v.It can be seen that the method is uniformly convergent and gives accurate results.It is clear from the figures and tables that the proposed method can handle the coupled Burgers' equations very accurately.
For finding an optimum value of the auxiliary term ℎ, the lowest point of the error of norm 2 of Function (40) is selected.The lowest point of  (ℎ), as h = 1.00000000005237, is gained by means of Maple mathematical software.By replacing the value of the auxiliary term in  (, , ℎ), the absolute error of the 11 th -order approximation of the proposed algorithm decreases remarkably.

Figure 1 .
Figure 1.Behavior of exact (Left) and approximate solutions (Right) of  for Test Problem 1.

Figure 2 .
Figure 2. Behavior of exact (Left) and approximate solutions (Right) of  for Test Problem 1.

Figure 1 .
Figure 1.Behavior of exact (Left) and approximate solutions (Right) of u for Test Problem 1.

Figure 1 .
Figure 1.Behavior of exact (Left) and approximate solutions (Right) of  for Test Problem 1.

Figure 2 .
Figure 2. Behavior of exact (Left) and approximate solutions (Right) of  for Test Problem 1.Figure 2. Behavior of exact (Left) and approximate solutions (Right) of v for Test Problem 1.

Figure 3 .
Figure 3.Comparison of absolute errors for different values of t for Test Problem 1.

Figure 3 .
Figure 3.Comparison of absolute errors for different values of  for Test Problem 1.

Figure 4 .
Figure 4. Approximation solution of u (Left) and approximation solution of v (Right) for different values of t.

Figure 4 .
Figure 4. Approximation solution of u (Left) and approximation solution of v (Right) for different values of t.

Figure 5 .
Figure 5. Behavior of exact (Left) and approximate solution (Right) for Test Problem 2.

Figure 6 .
Figure 6.Comparison of exact and approximate solution for Test Problem 2.

Figure 6 .
Figure 6.Comparison of exact and approximate solution for Test Problem 2.

Figure 6 .
Figure 6.Comparison of exact and approximate solution for Test Problem 2.

Figure 7 .
Figure 7.Comparison of the exact and approximate solution for  and  for  = 1 for Test Problem 3.

Figure 7 .
Figure 7.Comparison of the exact and approximate solution for u and v for t = 1 for Test Problem 3.

Table 1 .
Comparison of numerical results of u for different values of t and x for Test Problem 1.

Table 2 .
Comparison of numerical results for v for different values of t and x for Test Problem 1.

Table 3 .
Comparison of absolute errors for different values of t and x for Test Problem 1.

Table 4 .
Comparison of error norms different values of t and x for Test Problem 2.

Table 6 .
Comparison of L 2 and L ∞ error norms for different values of t for Test Problem 3.