Numerical similarity reductions of the (1+3)-dimensional Burgers equation

Applied Mathematics and Computation 217 (2011) 7455–7461

Contents lists available at ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

Numerical similarity reductions of the (1+3)-dimensional Burgers equation M.A. Christou a, C. Sophocleous b,⇑ a b

Department of Mathematics, University of Nicosia, CY 1678, Nicosia, Cyprus Department of Mathematics and Statistics, University of Cyprus, CY 1678, Nicosia, Cyprus

a r t i c l e

i n f o

Keywords: Burgers equation Similarity reductions Numerical solutions

a b s t r a c t We consider the (1+3)-dimensional Burgers equation ut = uxx + uyy + uzz + uux which has considerable interest in mathematical physics. Lie symmetries are used to reduce it to certain ordinary differential equations. We employ numerical methods to solve a number of these ordinary differential equations. Ó 2011 Elsevier Inc. All rights reserved.

1. Introduction We consider the (1+3)-dimensional Burgers equation

ut ¼ uxx þ uyy þ uzz þ uux :

ð1Þ

The one-dimensional case is the simplest and the best known equation in wave theory [18]. Furthermore it has applications in gas dynamics [6,14] and in plasma dynamics [22]. For further applications and exact solutions see for example Refs. [1,5,8,9,16,23,24]. Because of the wide range of applications of Burgers equation, several studies have been made to generalize it to higher dimensions. See for example [10] and the references therein for further applications of Eq. (1). Probably the most useful point transformations of partial differential equations are those which form a continuous (Lie) group of transformations, each member of which leaves an equation invariant. Symmetries of this equation are then revealed, perhaps suggesting links with equations studied in a different context, perhaps enabling new solutions to be found directly or via similarity reductions. Similarity reductions reduce the number of independent variables in a partial differential equation. In particular the reduction of a partial differential equation with respect to an r-dimensional (solvable) subalgebra of its Lie symmetry algebra leads to a reduction of the number of independent variables by r. Here the reduction of (1) can be done in three different ways using: (i) three-dimensional Lie subalgebras; (ii) two-dimensional Lie subalgebras; (iii) one-dimensional Lie subalgebras. The present paper is a continuation of the recent work [7]. The main goal of reference [7] was to derive exact solutions for Eq. (1) by reducing it to ordinary differential equations. This was achieved with the employment of Lie symmetries admitted by (1). The reduction was carried out in three different ways using: (i) three-dimensional Lie subalgebras; (ii) twodimensional Lie subalgebras; (iii) one-dimensional Lie subalgebras. In (i) Eq.(1) was directly reduced to an ordinary differential equation. In (ii) Eq. (1) was reduced to a partial differential in two independent variables. The next step was to derive the Lie symmetries for the reduced equation that led to transformations which mapped them into ordinary differential equations. Finally in (iii) Eq. (1) was reduced to a partial differential in three independent variables. The study of the Lie symmetries of the reduced equation produced further reductions. Lie symmetries of the latter reduced equations led to similarity transformations that mapped them into ordinary differential equations.

⇑ Corresponding author. E-mail addresses: [email protected] (M.A. Christou), [email protected] (C. Sophocleous). 0096-3003/$ - see front matter Ó 2011 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2011.02.042

7456

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

In the present paper we consider certain reduced ordinary and partial differential equations that we could not solve analytically in [7]. These equations are solved numerically using the Newton–Kantorovich method for nonlinear boundary value problems, [4,12,13,20]. The method is applied on both uniform and nonuniform meshes. In the current work we present numerical solutions obtained from three different cases, but the method can be applied in most cases reported in [7]. A more extensive description of the numerical method is given in Section 2. The procedure for finding Lie symmetries for differential equations is well known and established in the last few decades. See, for example, the textbooks [2,3,11,15,19]. The Lie symmetry algebra L of Eq. 1 is seven-dimensional and is spanned by the following symmetry operators [10]

v 1 ¼ @ t ; v 2 ¼ @ x ; v 3 ¼ @ y ; v 4 ¼ @ z ; v 5 ¼ 2t@ t þ x@ x þ y@ y þ z@ z  u@ u ; v 6 ¼ t@ x  @ u ; v 7 ¼ y@ z  z@ y :

ð2Þ

It is known that the reduction of a partial differential equation with respect to an r-dimensional (solvable) subalgebra of its Lie symmetry algebra leads to a reduction of the number of independent variables by r. Therefore in order to derive all possible reductions we need to classify the low-dimensional (one-, two- and three-dimensional) subalgebras of the maximal symmetry algebra L of Burgers equation 1. This goal is achieved using Ovsiannikov’s method of classification of subalgebras [21]. A brief description of this method and how the desired results are obtained can be found in [7] (see Theorem 1). 2. Description of the numerical approach Here we give a brief description of the applied numerical method [4,12,17,20]. In some cases we work directly on the reduced ordinary differential equations posed as two point boundary value problems using the Newton–Kantorovich method on a uniform mesh and in other cases we begin the investigation on the original partial differential equation and through the similarity reductions we are led to use the Newton–Kantorovich method on a nonuniform mesh. We consider the boundary value problem

  2 d F dF and a 6 x 6 b; ¼ g x ; F; dx2 dx

with b:c: FðaÞ ¼ a and FðbÞ ¼ b;

ð3Þ

where gðx; F; ddFxÞ is a nonlinear function of x, F and ddFx. The original problem, the (3D + 1) problem, is defined on the parallelogram

k1 6 x 6 k2 ; k3 6 y 6 k4 ; k5 6 z 6 k6 ðhere we consider : k1 ¼ k3 ¼ k5 ¼ 1; k2 ¼ k4 ¼ k6 ¼ 2Þ:

ð4Þ

We consider a uniform mesh on the original problem namely:

k2  k1 ; xi ¼ k1 þ ði  1Þhx ; Nx  1 k4  k3 ; yj ¼ k3 þ ðj  1Þhy ; hy ¼ Ny  1 k6  k5 hz ¼ ; zm ¼ k1 þ ðm  1Þhz ; Nz  1 Nx ¼ Ny ¼ Nz ¼ N is the number of grid points for each direction:

hx ¼

ð5Þ

When applied though the similarity transformations, the mesh becomes nonuniform in 1D. Therefore we need to define

hi ¼ xiþ1  xi ;

F 0i ¼

F iþ1  F i1

;

xiþ1  xi1 2½F iþ1 ðxi1xi Þ þ F i1 ðxi  xiþ1 Þ þ F i ðxiþ1  xi1 Þ 00 Fi ¼ ðxi1  xi Þðxi1  xiþ1 Þðxi  xiþ1 Þ

ð6Þ

and then the discrete form of our boundary value problem becomes

  2½F iþ1 ðxi1xi Þ þ F i1 ðxi  xiþ1 Þ þ F i ðxiþ1  xi1 Þ F iþ1  F i1 with F 1 ¼ a; F N ¼ b: ¼ g xi ; F i ; ðxi1  xi Þðxi1  xiþ1 Þðxi  xiþ1 Þ xiþ1  xi1

ð7Þ

In terms of the correction, D(x), the correction to the jth iteration of the Newton method for the problem is recast to

Dxx  g F 0 Dx  g F D ¼ g  ðF 00 Þj

and F jþ1 ðxÞ  F j ðxÞ þ DðxÞ:

ð8Þ

Now, if we apply the nonuniform mesh, we can observe that the approximation is not of order O(h2) as in the uniform grid but rather O(h). For example, if we consider the first derivative in a uniform mesh, we have 0

F0 ¼

2

3

0

2

3

00 000 00 000 F iþ1  F i1 F þ hF þ h2! F þ h3! F þ     ðF  hF þ h2! F  h3! F þ   Þ 2 ¼ F 0 þ Oðh Þ: ¼ 2h 2h

ð9Þ

7457

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

On the nonuniform mesh we consider the first derivative as

F iþ1  F i1

F0 ¼

xiþ1  xi1

¼

F iþ1  F i1

xiþ1  xi þ xi  xi1

¼

F iþ1  F i1 hi þ hi1

and we have 0 F iþ1  F i1 F þ hi F þ ¼ hi þ hi1

h2i 2!

h3

F 00 þ 3!i F 000 þ     ðF  hi1 F 0 þ hi þ hi1

h2i1 2!

F 00 

h3i1 3!

F 000 þ   Þ

¼ F0 þ

ðhi  hi1 Þ 00 F þ    ¼ F 0 þ OðhÞ: 2! ð10Þ

This obstacle does not interfere in the approximation or in the convergence of the algorithm because hi is of the order of 104, which is practically very small. This is because, if for example we consider in the original 3D problem 11 points for each direction, then in 1D this is translated to 113 = 1331 points and if we consider something like 31 points in 3D then we have 313 = 29791 points. This results the very small order of hi. In the case for which we treat directly the ordinary differential equation posed as a two-point boundary value problem on a uniform mesh we use (7) to obtain

2½hF iþ1 þ hF i1  2hF i  2h

3

    F iþ1  F i1 F iþ1 þ F i1  2F i F iþ1  F i1 ¼ g xi ; F i ; ¼ g x ; F ; ) with F 1 ¼ a; i i 2 2h 2h 2h

F N ¼ b: ð11Þ

3. Numerical solutions Recently [7] we have classified the list of similarity reductions that map (1) into ordinary differential equations. Here we consider certain reduced ordinary differential equations that cannot be solved analytically. We adopt the method described in the previous section to solve them numerically. As a first example, we consider the two-dimensional subalgebras of algebra 2 hv5 + a7v7, v1i. In order to find the corresponding similarity reduction we need to take two subcases: a7 = 0 and a7 – 0. Here we take a7 = 0. See Ref [7] for the subcase a7 – 0. The desired similarity reduction is obtained by solving the system

dt dx dy dz du ¼ ¼ ¼ ¼ ; 1 0 0 0 0

dt dx dy dz du ¼ ¼ ¼ ¼ : 2t x y z u

The three independent solutions of the system are I1 = ux, I2 ¼ yx and I3 ¼ xz. Hence the similarity reduction is

u¼

1 /ðg; nÞ; x

y x

g¼ ; n¼

z x

and it reduces (1) to

ð1 þ g2 Þ/gg þ 2gn/gn þ ð1 þ n2 Þ/nn þ 4g/g þ 4n/n  g//g  n//n  /2 þ 2/ ¼ 0: This partial differential equation admits the Lie symmetry hg@ n  n@ gi which, in turn, leads to the similarity transformation

/ ¼ FðxÞ;

x ¼ g2 þ n2 ¼

y2 þ z2 ; x2

which is obtained by solving the equation

dn

g

¼

dg d/ : ¼ 0 n

This mapping reduces the above partial differential equation to the ordinary differential equation

4ðx2 þ xÞF 00 þ ð10x  2xF þ 4ÞF 0  F 2 þ 2F ¼ 0:

ð12Þ

Here we consider the two point boundary value problem which consists of the above equation and the boundary conditions

Fðx1 Þ ¼ a1 ;

Fðx2 Þ ¼ b1 :

ð13Þ

To investigate the boundary value problem (12) and (13) we use the Newton–Kantorovich method on a nonuniform mesh as was described in the previous section. We used many cases with different pairs of boundary values and in all cases the method converged rapidly after using about ten iterations. The algorithm was blocked when the rate of convergence reached the value  = 1012. In Fig. 1 we demonstrate eight different solutions obtained using the same number of different pairs of boundary values. We solve and evaluate the function F(x) and then return and obtain u(x, y, z). Here we present sections of the solutions obtained by fixing the value for the variable x, say x = 1, and plot u(1, y, z) with respect to the other two variables y and z.

7458

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

-10

-6

-20

-10

-30

-14

1

1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2 1

-0.06

-2.0

-0.12

-2.6 1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2 1

1.1

4

1.04

3

1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2 1

2 1

10 8 6 4 1

1.2 1.4 1.6 1.8

1.8 1.6 1.4 1.2

2

5

1.16

1

2 1

2

0

-1.4

1

1.2 1.4 1.6 1.8

1.8 1.6 1.4 1.2

1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2

90 60 30

1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2 1

1.2 1.4 1.6 1.8

2 1

1.8 1.6 1.4 1.2

2

Fig. 1. 2-dimensional solutions for x = 1 obtained using different boundary conditions.

As a second example we consider the two-dimensional subalgebras of algebra 2 hv7 + a1v1, v2 + b1v1i. The solution of the corresponding system provides the similarity reduction

u ¼ /ðg; nÞ;

z

!

; g ¼ y2 þ z2 ; n ¼ t  b1 x  a1 sin1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 y þz

that maps (1) into 2

4g2 /gg þ ðb1 g þ a21 Þ/nn þ 4g/g  g/n  b1 g//n ¼ 0; where we take b1 – 0 since otherwise it is linear. If a1 – 0, this equation admits the Lie symmetry h@ ni which leads to the trivial solution u = c1ln (y2 + z2) + c2. If a1 = 0, it admits two-dimensional symmetry algebra with optimal system h@ ni, hb1n@ n + 2b1g@ g  (b1/ + 1)@ /i. The first component gives the trivial solution as before, while the second leads to the similarity transformation

/¼

  1 FðxÞ 1 ; b1 n

x¼

pffiffiffi

g

n

;

that reduces the above equation into the ordinary differential equation 2

2

2

ðb1 x3 þ xÞF 00 þ ð1 þ 4b1 x2 þ x2 FÞF 0 þ 2b1 xF þ xF 2 ¼ 0:

ð14Þ

7459

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

We note that (14) can be integrated once to give

1 2 2 ðb1 x3 þ xÞF 0 þ x2 F 2 þ b1 x2 F ¼ c1 : 2 An analytic solution exists for c1 = 0 [7]. Here we solve the boundary value problem which consists of the ordinary differential Eq. (14) and the boundary conditions

Fðx1 Þ ¼ a2 ;

Fðx2 Þ ¼ b2 :

ð15Þ

One needs to observe that for the given parameters the solution behaves like a traveling-wave solution with respect to the n variable. In Fig. 2 we present the solution obtained at four different times, for example t1 = 0, t2 = 0.5, t3 = 1 and t4 = 1.5. In this figure we plot the numerical solution of /(g, n) with respect to g and n. A better understanding of this graph is shown in Fig. 3 in which we present F(x) evaluated at t1, t2, t3 and t4. The initial solution is the first curve on the right-hand side of the figure and we see the shift to the left with the last solution obtained at t4 = 1.5. In the final example we consider a case from the optimal system of three-dimensional subalgebras of algebra 2. Here the similarity reduction maps the partial differential Eq. (1) directly into an ordinary differential equation. In particular we employ the element hv5, v3 + b2v2, v1i. We obtain the similarity reduction

u¼

1 FðxÞ; z

x¼

x  b2 y ; z

that transforms (1) into 2

ð1 þ b2 þ x2 ÞF 00 þ 4xF 0 þ 2F þ FF 0 ¼ 0: When we integrate once, we obtain

1 2 ð1 þ b2 þ x2 ÞF 0 þ 2xF þ F 2 ¼ k1 : 2 If k1 = 0, we derive the steady-state similarity solution, 2

4ð1 þ b2 Þ3=2 z : u ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffi h i 2 2 2 3=2 2y ffiffiffiffiffiffiffiffi þ 4k ð1 þ b Þ 1 þ b2 ðx  b2 yÞz þ ð1 þ b2 Þz2 þ ðx  b2 yÞ2 tan1 pxb 2 2 2

ð16Þ

1þb2 z

Now we consider the following boundary value problem evaluated on a uniform mesh 2

ð1 þ b2 þ x2 ÞF 00 þ 4xF 0 þ 2F þ FF 0 ¼ 0;

with FðaÞ ¼ a1 ;

FðbÞ ¼ a2 ;

a 6 x 6 b:

ð17Þ

In Fig. 4 we present the numerical solutions of (17) evaluated for different boundary conditions while in Fig. 5 we present solutions for the same problem calculated using an initial condition. In Fig. 5 we plot the numerical result and its initial

-1 -1.1 -1.2 -1.3 -1.4 -1.5

-1 -1.1 -1.2 -1.3 -1.4 -1.5

32 35 38 41 44 47 50

-3.25 -3.5 -3.75

-1 -1.1 -1.2 -1.3 -1.4 -1.5

32 35 38 41 44 47 50

-3.25

-3

-2.75

-1 -1.1 -1.2 -1.3 -1.4 -1.5

32 35 38 41 44 47 50

-4.25 -4.5 -4.75

32 35 38 41 44 47 50

Fig. 2. 2-dimensional solutions obtained using four different times.

-4.3

-4

-3.7

7460

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

2

1.5

1

0.5

0 -3

-2.8

-2.6

-2.4

-2.2

-2

-1.8

-1.6

-1.4

-1.2

-1

Fig. 3. 1-dimensional solutions obtained using four different times.

50 40 30 20 10 0 -10 -20 -30 -10

-7.5

-5.0

-2.5

0

2.5

5.0

7.5

10

Fig. 4. Numerical solutions of (17) obtained using different boundary conditions.

condition for four different cases. For example in the panel top-right we began with an initial condition the function sinx  cosx while in the panel placed bottom-left with the function 2

3

4ð1 þ b2 Þ2 qffiffiffiffiffiffiffiffiffiffiffiffiffiffi : 2 2 x þ 4k ð1 þ b2 Þ32 Þ 1 þ b2 x þ ð1 þ b2 þ x2 Þðarctan pffiffiffiffiffiffiffiffi 2 2 2 1þb2

We observe that the numerical method in this case is converging to the initial condition which means it was a very good initial guess, but here we need to state that the last initial condition is also the exact solution for a special case when k1 = 0, see (1).

7461

M.A. Christou, C. Sophocleous / Applied Mathematics and Computation 217 (2011) 7455–7461

140

1.5

120

1

100

0.5

80

0

60

-0.5

40

-1

20 0 -10

-5

0

5

10

18

0

2

4

6

8

10

1.2

16

1

14 12

0.8

10

0.6

8 6

0.4

4

0.2

2 0

-1.5

0

1

2

3

4

5

6

7

8

9

10

0 -10

-5

0

5

10

Fig. 5. Numerical solutions of (17) obtained using a given initial condition.

The same algorithm with some modifications can be applied to solve the remainder of the ordinary differential equations if we impose them as two-point boundary value problems. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24]

G.W. Bluman, J.D. Cole, The general similarity solution of the heat equation, J. Math. Mech. 18 (1968) 1025–1042. 69. G.W. Bluman, S. Kumei, Symmetries and Differential Equations, Springer, New York, 1989. G.W. Bluman, A.F. Cheviakov, S. Anco, Applications of Symmetry Methods to Partial Differential Equations, Springer, New York, 2009. J.P. Boyd, Chebyshev and Fourier Spectral Methods, second ed., Dover, New York, 2001. B.E. Clothier, J.H. Knight, I. White, Burgers’ equation: application to field constant-flux infiltration, Soil Sci. 132 (1981) 255–261. J.D. Cole, On a quasilinear parabolic equation occurring in aerodynamics, J. Appl. Math. 9 (1951) 225–236. M.A. Christou, N.M. Ivanova, C. Sophocleous, Similarity reductions of the (1+3)-dimensional Burgers equation, Appl. Math. Comput. 210 (2009) 87–99. E. Demetriou, M.A. Christou, C. Sophocleous, On the classification of similarity solutions of a two-dimensional diffusion–advection equation, Appl. Math. Comput. 187 (2007) 1333–1350. E. Demetriou, N.M. Ivanova, C. Sophocleous, Group analysis of (2 + 1)- and (3 + 1)-dimensional diffusion–convection equations, J. Math. Anal. Appl. 348 (2008) 55–65. M.P. Edwards, P. Broadbridge, Exceptional symmetry reductions of Burgers’ equation in two and three spatial dimensions, Z. Angew. Math. Phys. 46 (1995) 595–622. W.I. Fushchich, W.M. Shtelen, N.I. Serov, Symmetry Analysis and Exact Solutions of Equations of Nonlinear Mathematical Physics, Kluwer, Dordrecht, 1993. W.B. Gragg, R.A. Tapia, Optimal Error Bounds for the Newton–Kantorovich Theorem, SIAM J. Numer. Anal. 11 (1974) 10–13. S.N. Ha, A Nonlinear Shooting Method for Two-Point Boundary Value Problems, Comput. Math. Appl. 42 (2001) 1411–1420. E. Hopf, The partial differential equation ut + uux = luxx, Commun. Pure Appl. Math. 3 (1950) 201–230. N.H. Ibragimov, Elementary Lie Group Analysis and Ordinary Differential Equations, Wiley, New York, 1999. N.M. Ivanova, Exact solutions of diffusion–convection equations, Dyn. PDEs 5 (2008) 139–171. A. Klute, A numerical method for solving the flow equation for water in unsaturated materials, Soil Sci. 73 (1952) 105–116. M.J. Lighthill, Viscosity effects in sound waves of finite amplitude, in: G.K. Batchelor, R.M. Davis (Eds.), Surveys in Mechanics, Cambridge University Press, Cambridge, 1956, pp. 250–351. P. Olver, Applications of Lie Groups to Differential Equations, Springer-Verlag, New York, 1986. J.M. Ortega, The Newton–Kantorovich Theorem, Am. Math. Mon. 75 (1968) 658–660. L.V. Ovsiannikov, Group Analysis of Differential Equations, Academic Press, New York, 1982. P. Penel, C.M. Brauner, Identification of parameters in a nonlinear self-consistent system including a Burgers’ equation, J. Math. Anal. Appl. 45 (1974) 654–681. L. Rajaee, H. Eshraghi, R.O. Popovych, Multi-dimensional quasi-simple waves in weakly dissipative flows, Phys. D 237 (2008) 405–419. L.A. Richards, Capillary conduction of liquids through porous mediums, Physics 1 (1931) 318–333.