Wavelet-Based Homotopy Analysis Method for Nonlinear Matrix System and Its Application in Burgers Equation

Abstract

To generalize the homotopy analysis method (HAM) to multidegree-of-freedom nonlinear system, the adaptive precise integration method (APIM) is introduced into the HAM, with which the almost exact value of the exponential matrix can be obtained. Combining the interval interpolation wavelet collocation method, HAM-based APIM can be employed to solve the nonlinear PDEs. As an example, Burgers equation is spatially discretized by the interval quasi-Shannon wavelet collocation method and solved by the proposed method to illustrate the effectiveness and great potential of the homotopy analysis method in nonlinear problems.

1. Introduction

The homotopy analysis method (HAM) is developed in 1992 by Liao in [1–5]. This method has been successfully applied in various nonlinear problems in science and engineering in [6–11] and references therein. Based on homotopy of topology, the validity of the HAM is independent of whether or not there exist small parameters in the considered equation. Therefore, the HAM can overcome the foregoing restrictions and limitations of perturbation techniques [12, 13]. In fact, He [14, 15] proposed another asymptotic numerical method for solving nonlinear problems named homotopy perturbation method (HPM). The better improvement is adding an auxiliary parameter into the homotopy equation, which is helpful to eliminate the secular term in the perturbation solution. This can improve the rate of convergence greatly. Unlike analytical perturbation methods, HPM does not depend on small parameter which is difficult to find. Most researchers think of HPM as a special case of HAM, which is much simpler than classical HAM and has been widely used in solving various nonlinear problems [16]. However, the HAM applies mainly in single-degree-of-freedom systems; as for multidegree-of-freedom nonlinear differential equations, due to the correlation of the equations, the solving process is very tedious, and strong calculating technique is needed.

In an attempt to overcome the previous difficulties, the precise integration method (PIM) proposed by Zhong [17] in 1994 is employed, because it is effective in computation precision and stability for the solutions of linear differential equations [18, 19]. However, like other numerical methods, PIM has its own limitation. When it is applied to the problems of time-variant or nonlinear time integration, PIM usually expands the nonlinear term in Taylor series firstly, which is then solved by using an iteration method. In order to obtain high precision, a very small time step should be taken for PIM, which results in huge computational complexity and greatly restricts the application of PIM in engineering.

Burgers’ model of turbulence is a crucial fluid dynamic model, and the study of this model and theory of shock waves has been considered by many authors both for conceptual understanding of a class of physical flows and for testing various numerical methods. In this work, the Burgers equation is firstly transformed into a set of nonlinear matrix differential equations using wavelet collocation method [20–22]. Since the definition domain of wavelet transformation is an infinite interval, which can deduce the boundary effect, this can be eliminated by introducing an interval wavelet [23] instead of the common wavelet function.

Taking full advantage of the HAM and PIM, a novel approximate method for multidegree-of-freedom nonlinear matrix differential equation based on interval wavelet is presented in this paper.

2. The HAM Based on the PIM

2.1. Basic Idea of the Homotopy Analysis Method

Taking as an example, we consider the following nonlinear algebraic equation: 𝑓(𝑥)=0. (1) First of all, we construct such a homotopy 𝐻𝑥;𝑞=1−𝑞𝑓(𝑥)−𝑓𝑥₀+𝑞⋅𝑓(𝑥), (2) where 𝑥₀ is an initial guess of 𝑥 and 𝑞∈[0,1] is called embedding parameter. Obviously, at 𝑞=0 and 𝑞=1, one has 𝐻(𝑥;0)=𝑓(𝑥)−𝑓𝑥₀,𝐻(𝑥;1)=𝑓(𝑥), (3) respectively. Thus, as 𝑞 increases from 0 to 1, 𝐻(𝑥;𝑞) varies continuously from 𝑓(𝑥)−𝑓(𝑥₀) to 𝑓(𝑥). Such kind of continuous variation is called deformation in topology [20]. Now, enforcing 𝐻(𝑥;𝑞)=0, that is, 1−𝑞𝑓(𝑥)−𝑓𝑥₀+𝑞⋅𝑓(𝑥)=0, (4) we have now a family of algebraic equations. Obviously, the solution of the previous family of algebraic equations is dependent upon the embedding parameter 𝑞. So, the family of equations can be rewritten as 1−𝑞𝑓𝜙𝑥;𝑞−𝑓𝑥₀+𝑞⋅𝑓𝜙𝑥;𝑞=0. (5)

As the embedding parameter 𝑞 increases from 0 to 1, 𝜙(𝑥;𝑞) varies from the initial guess 𝑥₀ to the solution 𝑥 of 𝑓(𝑥)=0. We call the family of equations like (5) the zeroth-order deformation equation

Because 𝜙(𝑥;𝑞) is now a function of the embedding parameter 𝑞, we can expand it into Maclaurin series 𝜙𝑥;𝑞=𝑥₀+^+∞_𝑚=1𝑥_𝑚⋅𝑞^𝑚, (6) where 𝜙(𝑥;0)=𝑥₀ is employed, and 𝑥_𝑚=1⋅𝜕𝑚!^𝑚𝜙𝑥;𝑞𝜕𝑞^𝑚_𝑞=0. (7) Here, the series (6) is called homotopy series, and 𝑥_𝑚 is called the 𝑘th-order homotopy derivative of 𝜙(𝑥;𝑞). If the convergence radius of the previous Maclaurin series is no less than 1, using the relationship 𝜙(𝑥;1)=𝑥, one has 𝑥=𝑥₀+^+∞_𝑚=1𝑥_𝑚. (8)

According to the fundamental theorem of calculus about Taylor series, the coefficient 𝑥_𝑚 of the homotopy series (6) is unique.

Substituting (6) into (5) and rearranging based on powers of 𝑞-terms, one has 𝑥₁𝑓^′𝑥₀+𝑓𝑥₀𝑥=0,₂𝑓^′𝑥₀+12𝑥²₁𝑓^′′𝑥₀=0. (9) Whose solutions are 𝑥₁𝑓𝑥=−⁰𝑓^′𝑥₀,𝑥₂𝑥=−²¹𝑓^′′𝑥₀2𝑓^′𝑥₀𝑓=−²𝑥₀𝑓^′′𝑥₀2𝑓^′𝑥₀³. (10)(11)

In this way, one obtains 𝑥_𝑚 one by one in the order 𝑚=1,2,3,…. Then one has 1st-order homotopy-series approximation 𝑥≈𝑥₀+𝑥₁=𝑥₀−𝑓𝑥₀𝑓^′𝑥₀, (12) 2nd-order homotopy-series approximation 𝑥≈𝑥₀+𝑥₁+𝑥₂=𝑥₀−𝑓𝑥₀𝑓^′𝑥₀−𝑓²𝑥₀𝑓^′′𝑥₀2𝑓^′𝑥₀³. (13)

2.2. The Precise Integration Method for Nonlinear Matrix Differential Equations

Taking as an example, we consider the following matrix differential equations: 𝑑𝑣(),𝑑𝑡=𝐻𝑣+𝑓𝑣,𝑡 (14) where 𝑣 is an 𝑛-dimensional unknown vector, 𝐻 and 𝑓(𝑣,𝑡) are 𝑛×𝑛 matrix, and 𝑡 is a time variable. The general solution of (14) is 𝑣(𝑡)=𝑒^{𝐻(𝑡−𝑡𝑘)}𝑣𝑡_𝑘+𝑒^𝐻𝑡𝑡_𝑡𝑘𝑒^−𝐻𝑡⋅𝑓(𝑣,𝑡)𝑑𝑡. (15) To solve the integral part of (15), the linear interpolation approximation is employed 𝑓(𝑣,𝑡)=𝑟₀+𝑟₁𝑡−𝑡_𝑘, (16) where 𝑟₀=𝑓(𝑣(𝑡_𝑘),𝑡_𝑘), 𝑟₁=𝑟₁=(𝜕𝑓/𝜕𝑡)|_{𝑡=𝑡𝑘}+∑^2𝑛_𝑖=1[(𝜕𝑓/𝜕𝑣_𝑖)(𝐻𝑣(𝑡_𝑘)+𝑓(𝑣(𝑡_𝑘),𝑡_𝑘))_𝑖], and 𝑣_𝑖 is the 𝑖th element of vector 𝑣. Based on the previous derivations, the recurrence formula of (14) can be obtained as 𝑣𝑡_𝑘+1=𝑒^{𝐻(𝑡𝑘+1−𝑡𝑘)}𝑣𝑡_𝑘+𝐻⁻¹𝑟₀+𝐻⁻¹𝑟₁−𝐻⁻¹𝑟₀+𝐻⁻¹𝑟₁+𝑡−𝑡_𝑘𝑟₁. (17)

Supposing the time step 𝜏=𝑡_𝑘+1−𝑡_𝑘 is constant, the matrix exponential function can be calculated accurately by using PIM. The calculation of the exponential function is based on the identity 𝑒^𝐻⋅𝜏=𝑒^{𝐻⋅(𝜏/𝑚)𝑚}, (18) where 𝑚 is an integer. It is suggested to select 𝑚=2^𝑁, such as 𝑁=20, 𝑚=1048576.

Because 𝜏 is a small time interval, Δ𝑡=𝜏/𝑚 is an extremely small time interval. Hence for the interval, the truncated Taylor expansion is applied with high precision 𝑇=𝑒^{𝐻⋅Δ𝑡}=𝐼+𝑇𝑎, (19) where 𝑇𝑎=𝐻⋅Δ𝑡+(𝐻⋅Δ𝑡)²/2+(𝐻⋅Δ𝑡)³/3!+(𝐻⋅Δ𝑡)⁴/4!

Because Δ𝑡 is very small, the first five-term series expansion should be enough.

For computing the matrix 𝑇, (17) could be factored as 𝑇=[𝐼+𝑇𝑎]^2𝑁=[𝐼+2⋅𝑇𝑎+𝑇𝑎⋅𝑇𝑎]^2𝑁−1. (20) Such factorization should be iterated 𝑁 times. Noted that it is equivalent to do the following instruction: for(iter=0;iter!=𝑁;++iter)𝑇𝑎=2⋅𝑇𝑎+𝑇𝑎⋅𝑇𝑎, (21) after the execution of the above instruction, 𝑇=𝐼+𝑇𝑎.

Because the nonlinear term 𝑓(𝑣,𝑡) is expanded in Taylor series and the PIM uses an iterative method, the time step 𝜏 must be very small.

3. Coupling Technique of the HAM and the PIM for Burgers Equation Based on Interval Wavelet

3.1. Spatial Discretization of the Burgers Equation

Consider that the one-dimensional Burgers equation has the form [21–23] 𝜕𝑢𝜕𝑡+𝑢𝜕𝑢=1𝜕𝑥𝜕Re²𝑢𝜕²𝑥,𝑥∈[0,2],𝑡≥0, (22) subject to the initial and boundary conditions 𝑢(𝑥,0)=sin(𝜋𝑥),𝑢(0,𝑡)=𝑢(2,𝑡)=0, (23) where 𝑡 represents the time, while Re denotes the Reynolds number.

According to the basic idea of wavelet collocation method, the solution 𝑢(𝑥,𝑡) could be approximately represented as 𝑢(𝑥,𝑡)≈𝑢_𝑗(𝑥,𝑡)=^2𝑗_𝑛=0𝑢𝑥_𝑛,𝑡𝜑_𝑗𝑥−𝑥_𝑛, (24) where 𝑥_𝑛(𝑛=0,1,2,…,2^𝑗) are uniform discrete points and 𝜑(𝑥) is the scaling function. Here we take the quasi-Shannon wavelet [24, 25] as the scaling function 𝜑(𝑥)=sin(𝜋𝑥)−𝑥𝜋𝑥exp²2𝜎²,𝜎>0, (25) where 𝜎 is the window size parameter.

In order to take the quasi-Shannon scaling function as the basis function, 𝑢(𝑥,𝑡) is discretized in its definition domain 𝑥_𝑛=2^1−𝑗⋅𝑛,𝑗∈𝑍,𝑛=0,1,2,…,2^𝑗. (26) And thereby the basis function is obtained 𝜑_𝑗𝑥−𝑥_𝑛=sin2^𝑗−1𝜋𝑥−𝑥_𝑛2^𝑗−1𝜋𝑥−𝑥_𝑛−2exp^2𝑗−3𝑥−𝑥_𝑛²𝑟², (27) where 𝑟=2^𝑗−1⋅𝜎.

It is easy to see that 𝜑_𝑗(𝑥−𝑥_𝑛) has the interpolation property 𝜑_𝑗𝑥_𝑘−𝑥_𝑛=𝛿_𝑘,𝑛=1𝑘=𝑛0𝑘≠𝑛. (28)

Using its interpolation property, the first derivative and second derivative of 𝜑_𝑗(𝑥−𝑥_𝑛) can be obtained as follows: 𝜑^′_𝑗𝑥−𝑥_𝑛=20𝑘=𝑛_𝑗−1cos[𝜋(𝑘−𝑛)]exp−(𝑘−𝑛)²/2𝑟²()𝜑𝑘−𝑛𝑘≠𝑛,_𝑗′′𝑥−𝑥_𝑛=−3+𝜋²𝑟²⋅2^2𝑗−23𝑟²𝑘=𝑛−2^2𝑗−2⋅1(𝑘−𝑛)²+1𝑟²[𝜋(−()⋅cos𝑘−𝑛)]⋅exp𝑘−𝑛²2𝑟²𝑘≠𝑛. (29)

Substituting (4) into (1), a series of nonlinear ODEs can be obtained 𝑑𝑢𝑥_𝑘,𝑡=𝑑𝑡_2𝑗𝑛=0𝑢𝑥_𝑛1,𝑡𝜑Re^𝑗′′𝑥_𝑘−𝑥_𝑛−𝑢𝑥_𝑘,𝑡𝜑^′_𝑗𝑥_𝑘−𝑥_𝑛. (30)

For simplicity, we can write (10) as matrix differential equation 𝑑𝑉(𝑡)𝑑𝑡=𝑀₀⋅𝑉(𝑡)+𝑀₁(𝑉(𝑡))⋅𝑉(𝑡), (31) where𝑉(𝑡)=𝑢𝑥₀,𝑡,𝑢𝑥₁,𝑡,…,𝑢𝑥₂^𝑗,𝑡^𝑇,𝑀₀=1𝜑Re^𝑗′′𝑥₀−𝑥₀𝜑_𝑗^′′𝑥₀−𝑥₁⋯𝜑_𝑗^′′𝑥₀−𝑥₂^𝑗𝜑_𝑗^′′𝑥₁−𝑥₀𝜑_𝑗^′′𝑥₁−𝑥₁⋯𝜑_𝑗^′′𝑥₁−𝑥₂^𝑗𝜑⋮⋮⋯⋮_𝑗′′𝑥₂^𝑗−𝑥₀𝜑_𝑗^′′𝑥₂^𝑗−𝑥₁⋯𝜑_𝑗^′′𝑥₂^𝑗−𝑥₂^𝑗,𝑀₁(𝑉(𝑡⋅𝜑=−diag))_′𝑗𝑥₀−𝑥₀𝜑^′_𝑗𝑥₀−𝑥₁⋯𝜑^′_𝑗𝑥₀−𝑥₂^𝑗𝜑^′_𝑗𝑥₁−𝑥₀𝜑^′_𝑗𝑥₁−𝑥₁⋯𝜑^′_𝑗𝑥₁−𝑥₂^𝑗𝜑⋮⋮⋮_′𝑗𝑥₂^𝑗−𝑥₀𝜑^′_𝑗𝑥₂^𝑗−𝑥₁⋯𝜑^′_𝑗𝑥₂^𝑗−𝑥₂^𝑗. (32)

3.2. Construction of Quasi-Shannon Interval Wavelet

The Lagrange polynomials 𝑙¹_𝑗,𝑘, 𝑙²_𝑗,𝑘 are introduced 𝑙¹_𝑗,𝑘=^𝐿_{𝑖=0𝑖≠𝑘}𝑥−𝑥_𝑗,𝑖𝑥_𝑗,𝑘−𝑥_𝑗,𝑖,𝑙²_𝑗,𝑘=^2𝑗_{𝑖=2𝑗−𝐿𝑖≠𝑘}𝑥−𝑥_𝑗,𝑖𝑥_𝑗,𝑘−𝑥_𝑗,𝑖. (33) The following two weight coefficients are defined 𝑎_𝑛𝑘=𝑙¹_𝑗,𝑘𝑥_𝑗𝑛,𝑏_𝑛𝑘=𝑙²_𝑗,𝑘𝑥_𝑗𝑛. (34) Then the interpolation basic function can be expressed as 𝑤_𝑗,𝑘(𝑥)=𝑤_𝑗,𝑘(𝑥)+⁻¹_{𝑛=−𝑁+1}𝑎_𝑛𝑘⋅𝑤_𝑗,𝑛(𝑥)𝑤,𝑘=0,…,𝐿,_𝑗,𝑘(𝑥)=𝑤_𝑗,𝑘(𝑥),𝑘=𝐿+1,…,2^𝑗𝑤−𝐿−1,_𝑗,𝑘(𝑥)=𝑤_𝑗,𝑘(𝑥)+^{2𝑗+𝑁−1}_𝑛=2𝑗+1𝑏_𝑛𝑘⋅𝑤_𝑗,𝑛(𝑥),𝑘=2^𝑗−𝐿,…,2^𝑗, (35)(36)(37) where 𝐿 is the number of he external points and 𝑁 is the support domain of the wavelet function.

It is easy to know that the quasi-Shannon interval wavelet is a linear combination of the quasi-Shannon scaling function 𝑤_𝑗,𝑘(𝑥). Therefore, the quasi-Shannon interval wavelet function possesses all the properties of the quasi-Shannon scaling function.

3.3. The PIM for Matrix Differential Equations

We choose the time step as 𝜏, assuming that the nonlinear term in (31) is linear in the relatively short time segment (𝑡_𝑘,𝑡_𝑘+1), the nonlinear term 𝑀₁(𝑉(𝑡))⋅𝑉(𝑡) of (31) can be expressed in Taylor series, and the first two terms are taken. Therefore, (31) could be represented as 𝑑𝑉(𝑡)𝑑𝑡=𝑀₀𝑉(𝑡)+𝑟₀+𝑟₁𝑡−𝑡_𝑘, (38) where 𝑟₀=𝑀₁𝑡_𝑘𝑉𝑡_𝑘,𝑟₁=𝑀₂𝑡_𝑘+𝑀₁𝑡_𝑘𝑀₀+𝑀₁𝑡_𝑘𝑀₁𝑡_𝑘𝑉𝑡_𝑘,𝑀₂(𝑡)=𝑑𝑀₁(𝑡).𝑑𝑡 (39)

Based on the previous derivation, the iterative solution of (38) can be obtained 𝑉𝑡_𝑘+1=𝑒^𝑀0⋅𝜏𝑉𝑡_𝑘+𝑀₀⁻¹𝑟₀+𝑀₀⁻¹𝑟₁−𝑀₀⁻¹𝑟₀+𝑀₀⁻¹𝑟+𝑟₁𝜏, (40) where the matrix exponential function 𝑒^𝑀0⋅𝜏 can be calculated accurately by using PIM.

3.4. The Homotopy Analysis Method Based on the Precise Integration Method

There are various ways to construct a homotopy function. For (11), a linear homotopy can be constructed as 𝑞𝑑𝜙𝑥,𝑡;𝑞𝑑𝑡−𝑀₀−𝑀₁𝜙𝜙+×𝑥,𝑡;𝑞𝑥,𝑡;𝑞1−𝑞𝑑𝜙𝑥,𝑡;𝑞𝑑𝑡−𝑀₀𝜙𝑥,𝑡;𝑞−𝑀₁(𝐴)𝜙𝑥,𝑡;𝑞=0, (41) where 𝐴=𝑉(0) is the initial approximation of 𝑉(𝑡) and 𝑞∈[0,1] is the embedding parameter. As 𝑞 increases from 0 to 1, the solution 𝜙(𝑥,𝑡;𝑞) varies from the initial guess 𝐴 to the solution 𝑉(𝑡), and it can be expanded in Taylor series about the embedding parameter 𝜙𝑥,𝑡;𝑞=𝐴+^+∞_𝑚=1𝑉_𝑚(𝑡)⋅𝑞^𝑚, (42) where 𝑉_𝑚(𝑡)=1𝜕𝑚!^𝑚𝜙𝑥,𝑡;𝑞𝜕𝑞^𝑚 (43)

that the approximate solutions are generally not valid for large 𝑡. In order to improve the computation accuracy, we divide the time interval [0−𝑡] into a number of sections uniformly by setting the time step as 𝜏, and use the iterative method. In such a way, a series of time nodes can be expressed as 𝑡₀=0,𝑡₁=𝜏,…,𝑡_𝑘=𝑘⋅𝜏,…. (44)

Since the only requirement for 𝑞 is that from 𝑞=0 to 𝑞=1, 𝜙(𝑥,𝑡;𝑞) should change continuously from the initial guess 𝐴 to 𝑉(𝑡). Therefore, it is feasible to improve the computation accuracy by assuming that 𝑞=Δ𝑡𝜏,Δ𝑡∈[0,𝜏]. (45) Then 𝜙(𝑥,𝑡;𝑞) can be expanded in the form 𝜙𝑥,𝑡;𝑞=𝑉₀(𝑡)+𝑉₁(𝑡)⋅Δ𝑡𝜏+𝑉₂(𝑡)⋅Δ𝑡𝜏²+𝑉₃(𝑡)⋅Δ𝑡𝜏³+⋯. (46)

Substituting (46) into (41) and rearranging based on powers of 𝑞-terms, the following are obtained: 𝑉₀(𝑡)=𝐴,𝑑𝑉₀(𝑡)𝑑𝑡=𝐻⋅𝑉₀(𝑡)𝑉=𝐻⋅𝐴,₁(𝑡)=𝜏⋅𝑑𝑉₀(𝑡)𝑑𝑡=𝜏⋅𝐻⋅𝐴,𝑑𝑉₁(𝑡)𝑑𝑡=𝐻⋅𝑉₁(𝑡)₀=𝜏⋅𝐻²𝑉⋅𝐴,₂(𝑡)=𝜏2⋅𝑑𝑉₁(𝑡)=()𝑑𝑡𝜏𝐻²2⋅𝐴,𝑑𝑉₂(𝑡)𝑑𝑡=𝐻⋅𝑉₂(𝑡)+𝑠₀=𝜏²𝐻³2⋅𝐴+𝑠₀,𝑉₃(𝑡)=𝜏3⋅𝑑𝑉₂(𝑡)=()𝑑𝑡𝜏𝐻³𝜏3!⋅𝐴+3⋅𝑠₀,𝑑𝑉₃(𝑡)𝑑𝑡=𝐻⋅𝑉₃(𝑡)+𝑠₁=𝜏³𝐻⁴3!⋅𝐴+𝜏𝐻3⋅𝑠₀+𝑠₁,𝑉₄(𝑡)=𝜏4⋅𝑑𝑉₄(𝑡)=()𝑑𝑡𝜏𝐻⁴𝜏4!⋅𝐴+²𝐻3⋅4⋅𝑠₀+𝜏4⋅𝑠₁, (47) where 𝑠₀=𝑀₁(𝑉₁(𝑡))⋅𝐴, 𝑠₁=𝑀₁(𝑉₁(𝑡))⋅𝑉₁(𝑡)+𝑀₁(𝑉₂(𝑡))⋅𝑉₀(𝑡) and 𝐻=𝑀₀+𝑀₁(𝐴).

Substituting (47) into (46) and assuming 𝑞=1, one has 𝑉(𝑡)=𝑉₀(𝑡)+𝑉₁(𝑡)+𝑉₂(𝑡)+𝑉₃(𝑡)=()+⋯𝐼+𝜏𝐻+𝜏𝐻²𝐴+𝜏2!+⋯3+𝜏²𝐻+𝜏3⋅4³𝐻²𝑠3⋅4⋅5+⋯₀+𝜏4+𝜏²𝐻+𝜏4⋅5³𝐻²𝑠4⋅5⋅6+⋯₁𝐻+⋯=𝑇⋅𝐴+2⋅⁻¹³𝜏²(𝜏𝐻)³+()3!𝜏𝐻⁴𝑠4!+⋯₀𝐻+3!⋅⁻¹⁴𝜏³(𝜏𝐻)⁴+()4!𝜏𝐻⁵𝑠5!+⋯₁𝐻+⋯=𝑇⋅𝐴+2⋅⁻¹³𝜏²⋅𝑅⋅𝑠₀𝐻+3!⋅⁻¹⁴𝜏³⋅𝑄⋅𝑠₁+⋯, (48) where 𝑅=[𝑇−𝐼−𝜏𝐻−(𝜏𝐻)²/2], 𝑄=𝑅−(𝜏𝐻)³/6, 𝑇=𝑒^𝜏𝐻. The matrix exponential function 𝑇 can be calculated by using PIM.

For simplicity, only the first two terms are taken as the approximate solution. Accordingly, the approximate solution is obtained 𝑉(𝑡)≈𝑉₀(𝑡)+𝑉₁(𝑡)𝐻=𝑇⋅𝐴+2⋅⁻¹³𝜏²⋅𝑅⋅𝑠₀. (49)

4. Numerical Results and Discussion

The Burgers equation has analytical solution as follows:𝑢()=𝑥,𝑡4𝜋⋅Re^∞_𝑛=1−𝑛𝑛exp²𝜋²𝑡𝐼Re^𝑛Re()×22𝜋sin𝑛𝜋𝑥_∞𝑛=1−𝑛𝑛exp²𝜋²𝑡𝐼Re^𝑛Re()2𝜋cos𝑛𝜋𝑥+𝐼₀Re2𝜋⁻¹, (50)where 𝐼_𝑛(𝑥) is the 𝑛th order modified Bessel function of the first kind. Figure 1 shows the analytical solutions with Re = 100 at 𝑡=0, 0.4, 0.6, respectively. It is easy to see that the analytical solution evolves into a shock wave near 𝑥=1, and the gradient becomes larger and larger with the increasing of Re.

Figure 1 

Analytical solutions of the Burgers equation at different times (𝑡=0,0.4,0.6).

In the analysis of the numerical result, we choose Reynolds number Re=100, the scaling parameter 𝑗=8, the artificial parameter 𝑟=3, and the number of the external points 𝐿=4.

The comparison of the classical PIM and the HAM is shown in Figures 2 and 3 and Tables 1 and 2.

Table 1 

Difference between the classical PIM and the HAM when 𝜏=0.001 and 𝜏=0.005.

Table 2 

Difference between the classical PIM and the HAM when 𝜏=0.05 and 𝜏=0.1.

(a) 𝜏=0.05,  𝑡=0.6

(b) 𝜏=0.1,  𝑡=0.6

Figure 2 

Numerical results of the classical PIM and the HAM.

(a) 𝜏=0.005,  𝑡=0.6

(b) 𝜏=0.001,  𝑡=0.6

Figure 3 

Errors of the classical PIM and the HAM.

We can see clearly that when the time step 𝜏 is relatively large (𝜏=0.05 or 𝜏=0.1), the classical PIM fails actually it fails when 𝜏 is 0.04. However, the HAM has enough precision even when the time step 𝜏 increases to 0.1 and 𝑡 is more than 0.6.

As the time step 𝜏 becomes smaller, both the classical PIM and the HAM have higher precision; however, we can see clearly in Figure 3 that the HAM is still better than the classical PIM. Furthermore, from the comparison of (49) and (40), it can be seen that the computation cost of the proposed method is less than the classical PIM.

5. Conclusions

In this paper, it is shown that the proposed method can solve the Burgers equation successfully. The numerical comparisons presented demonstrate that the solutions obtained by the classical PIM are not valid for large time step; the method proposed in this paper can obtain higher accuracy than the classical PIM even if the time step is small. Although the method proposed aims at Burgers equation, it is a promising tool for other nonlinear problems.

Acknowledgment

This work is supported by the National Key Technologies R & D Program of China under Grant no. 2012BAD35B02 and the National Natural Science Foundation of China under Grant no. 41171337.