A Mathematical Approach Based on the Homotopy Analysis Method: Application to Solve the Nonlinear Harry-Dym (HD) Equation

In this paper, the homotopy analysis method (HAM) has been employed to obtain the approximate analytical solution of the nonlinear Harry-Dym (HD) equation, which is one of the most important soliton equations. Utilizing the HAM, thereby employing the initial approximation, variations of the 7th-order approximation of the Harry-Dym equation is obtained. It is found that effect of the nonzero auxiliary parameter on convergence rate of the series solution is undeniable. It is also shown that, to some extent, order of the fractional derivative plays a fundamental role in the prediction of convergence. The final results reported by the HAM have been compared with the exact solution as well as those obtained through the other methods.


Introduction
The nonlinear Harry-Dym equation is [1], with the initial approximation, ( ) where a and b are constants.Equation ( 1) is one of the most important soliton equations i.e., a soliton is

Homotopy Analysis Method (HAM)
The essence of the HAM [15]- [20] is based on the approximation of PDEs along with a convenient way to guarantee convergence of the series solution.This issue can be outlined by the following example.
Consider the following nonlinear differential equation, ( ) in which A is a nonlinear operator, x and t are independent variables and ( ) , u x t is an unknown function.On account of a continuous mapping from ( ) ( ) , ; , , .
which must be one of the solutions of Equation ( 3), as proved by Liao [15].Now, it is of particular interest to show that Equation ( 7) provides a relationship between the initial approximation and the exact solution so that: 1) Solution of Equation ( 4) is valid for 2) The deformation derivative exists for all values of m.
Differentiating the zeroth-order deformation Equation ( 4) m times with respect to p, then setting 0 p = and dividing them by !m , gives the so-called mth-order deformation equation as, ) It should be noted that ( ) , m u x t for 1 m ≥ is governed by the linear Equation (8) with the linear boundary conditions that come from the original problem, which can be easily solved using some software programs such as MAPLE or MATHEMATICA.
Then can HAM be expeditiously used in NPDEs?There have been a number of NPDEs which can be solved through the HAM.To address this subject, Nassar et al. [21] developed the HAM to derive an approximate solution of the Poisson-Boltzmann equation for semiconductor devices.The solutions of the fractional Swift Hohenberg equation can be found in the work of Vishal et al. [22].They could report influence of the real bifurcation parameter on probability density function and then concluded that the proper values of auxiliary and homotopy parameters are needed.However, their model was not capable of taking into account the two-dimensional domain for fractional Brownian motion.
Hetmaniok et al. [23]  integral equations of the second kind.They also used the auxiliary parameter to show influence of this parameter on the convergence rate.In a book, current advances of the HAM to solve strongly nonlinear problems were edited by Liao [24].Ren et al. [25] introduced an approximate analytical solution to the Gross-Pitaevskii equation while a nonlinear Schrödinger equation was utilized to simulate Bose-Einstein condensates trapped in a harmonic potential.Martin [26] proposed an analytical algorithm to solve an integral-differential equation of the transport theory for stationary case using HAM.The solution of nonlinear boundary value problems (NBVPs) based on the Chebyshev operational matrices was investigated by Shaban et al. [27].They applied the Tau method to convert a set of algebraic equations so that the solution can be obtained iteratively.
Motsa et al. [28] utilized the spectral-homotopy analysis method (SHAM) to solve the Darcy-Brinkman-Forchheimer equation and found the admissible range of auxiliary parameter for this problem.Hesameddini and Latifizadeh [29] obtained solutions for the first and second types of the Painlevé equations with 8, 10 and 12 iterations.Wang [30] developed an approximate analytical solution of relativistic Toda lattice system and concluded that using HAM may be very effective for solving differential difference equations (DDEs).In other words, some papers are devoted to compare the HAM with the other analytical methods e.g., optimal homotopy asymptotic method (OHAM), differential transformation method (DTM), homotopy perturbation method (HPM), general series expansion method, harmonic balance method (HBM) which can be found in Refs.

HAM Formulation for the Nonlinear Harry-Dym Equation
As mentioned above, the Equation (1) can be solved through the HAM.To this end, the linear operator takes the form of, with the property, in which d is an integral constant.Now the nonlinear operator can be expressed in terms of Equation (1) as, Using the above definition, with the assumption ( ) , 1 H x t = , the zeroth-order deformation can be given as, For 0 p = and 1 p = , it is respectively written as, Herein, the th-order deformation equation is given by, in which, ( ) ( ) ( ) ( ) It should be noted that the discretized form of above equation is illustrated in Appendix B. So, solution of the th-order deformation equation for 1 m ≥ can be expressed by, Now, according to Equation ( 18), the values of ( ) ( ) ( ) To continue this procedure for 3 m > , ( ) , m u x t can be achieved and the series solution is completely computed.

Comparison and Validation
To validate the present iterative approach, results are compared with the exact solution of the Harry-Dym equation provided by Mokhtari [36] which can be expressed as, ( ) ( ) where c is a constant.Mokhtari [36] also showed that the variational iteration method (VIM) is not able to obtain implicit solution of the Harry-Dym equation.
However, one should point out that Hereman et al. [37] proposed implicit solution of this equation via a direct method and also performed transformations between the Harry-Dym and KdV equations explicitly.In particular, according to the fractional Harry-Dym equation, where 0 1 α ≤ ≤ is a parameter describing order of the fractional derivative, predicted by the HAM is in good agreement with those provided by Mokhtari [36].However, there is a relative error between the results which can be calculated by, where ( ) , u x t and ( ) , u x t are the exact and approximate solutions of Equation (1), respectively.Due to 4 a = and 1 b = [36], it is noted that the maxi- mum error does not exceed 1.098%, as shown in Figure 2. The root-meansquare error (RMSE) between the exact solution and the present iterative  approach, i.e., ( ) ( ) ( ) , , is also 0.509%.Apart from the exact solution, other observations may be made from Figure 1; e.g., it is assumed that for 1 h = − , the series solution will be converged at the results presented in [38] [39], and in general, the HAM consists of the HPM and Adomian decomposition method (ADM).In fact, the two exhibit similar pattern for different values of the auxiliary parameter.Therefore, it is particularly important to account for the interaction between the auxiliary parameter and initial approximation before using the HAM even for the cases which NPDEs have been taken into account.This type of scheme has been presented by Turkyilmazoglu [43] as an illustrative comparison between the HAM and HPM.In addition, using 7thorder approximation of the HAM instead of ADM presented in Ref. [44] is sug-gested largely due to its less time-consuming and also is able to obtain valid values of the auxiliary parameter.However, it should be noted here that although we utilize 1 h = − in our solution, the convergence can be observed when the objective is to replace Equation ( 21) by Equation (1) for Another comparison is also performed here to illustrate the accuracy of present iterative approach when the effect of nonlinearity has been taken into account.In connection with this matter, Rawashdeh [40] utilized the fractional reduced differential transformation method (FRDTM) to solve the nonlinear Harry-Dym equation while the maximum error between the HAM and FRDTM for ( ) , u x t is 11.723%.In addition, He found that ( ) , u x t predicted by the FRDTM is on average 1.18 less than that predicted by exact solution and is greatly affected by the number of iterations performed.Based on this comparison, it is seen that Equation ( 21) without 1 α = would also fit into HAM, namely, Caputo sense [45], which is not correspondent to define the fractional order initial conditions.Herein, the Caputo's formula is defined as [46], It should be noted that Kumar and Singh [41] applied the Laplace transform on the both sides of Equation ( 21) and solved the governing equation through a homotopy perturbation transform method (HPTM), but their solution could not capture ( ) , u x t especially for the cases with large values of auxiliary parameter.There exists a discrepancy between the exact solution and those obtained through HPTM because of decomposing the nonlinear term in He's polynomial (i.e., Equation ( 17)).To clarify, the He's polynomial has been also used in some previous studies [47] [48].It is to be noted that Guo and Mei [47] studied timefractional Boussinesq-type equation by the VIM and Khan and Wu [48] investigated the advection problems through HPTM.

Universal Graphs for Convergence Region of the Presented Iterative Approach
The fact that Equation ( 18) is not convergent for all values of the auxiliary parameter, lifts some particular restrictions on choosing a proper value of , h h h ∈ (e.g., by trial and error [49]).Above all, making a good selection of h is important not only with the convergence, but also with the basis functions [16] [50].It should be noted that other perturbation methods require the need for small parameters, which may commonly diverge series solutions with slight changes (see Ref. [17]).Indeed, at large values of the auxiliary parameter, convergence may lead to only few terms of series solution.It is shown that based on the zerothorder deformation Equation ( 4), one may expect the flexibility and freedom to select L, despite the fact that the order of linear operator will be different from the original nonlinear problem (see Ref. [50]).As said in Ref. [49] "••• we are free to enhance the convergence region and convergence rate of a series solution via an appropriate choice of the auxiliary parameter even for -R.A. Van Gorder and K. Vajravelu [49] Due to Figure 1 that ( ) , u x t varies almost linearly versus x, Figure 3 shows the variation of ( ) , u x t for some different auxiliary parameters.It is seen that 1 h < − indicates the ascending behavior while 1 h > − suggests the descending behavior of the present iterative approach.Although the slopes of these lines are slightly different from each other, the lines become parallel by increasing or decreasing the auxiliary parameter in the neighbourhood of 1 h = − .It should be noted that this fact can result in almost linear variation of ( ) , u x t .In theory, at the mth-order of approximation, the square residual error can be defined as [51], It should be noted that as m ∆ decreases to zero more rapidly, the conver- gence rate for corresponding series solution would become faster [51].In this regard, Yabushita et al. [52] introduced an analytical approximation to projectile motion with the quadratic resistance law by the HAM.They derived governing equations using a power series which was accurate enough for large angle of projection.Therefore, with a minimum value of m ∆ , the optimal values of un- known auxiliary parameters were obtained.It should be noted that by increasing t, ( ) , u x t will be decreased slightly, as it is shown in Figure 4(a) and Figure 4(c).It is mainly due to the fact that this parameter plays an important role in higher-order governing equations (e.g., Equation (19c)) (see Refs. [22] [33] [53] [54]).It is to be noted that a similar conclusion for increasing t can also be drawn through a similar procedure.Odibat [55] reached a pragmatic approach for the convergence of HAM, and to date, no counterexample exits in the open literature for this conclusion, neither for NPDEs nor for PDEs.Hence, in such graphs we increase t from 1 to its a neighborhood which can exhibit this deviation for a better presentation.However, it is to be noted that this parameter is small in most of the analytical approximations (e.g., 0 2 t ≤ ≤ ) [18] [19] [49].
According to the higher-order equations (e.g., Equation (19c)), if t will be increased, the maximum value of ( ) , u x t will be increased for 1 1 h − ≤ ≤ (see Figure 4(b) and Figure 4(d)).It is noted that ( ) , u x t almost drops to zero when the auxiliary parameter reaches to zero as well.It is seen that neglecting the effect of t may lead to underestimation of ( ) , u x t versus h.However, if this parameter exceeds its 1 t , which can be named as a minimum value of t, ( ) , u x t will be fully negative.
It is noteworthy to mention that although many techniques (e.g., HPM, ADM, FRDTM, HPTM) can be useful for predicting analytical approximations of the nonlinear Harry-Dym Equation ( 1), the more accurate results can be obtained  through the HAM.On the other hand, since the HAM is employed to approximate NPDEs analytically, it is important to note that it enables us to construct a mapping of initial approximation to the exact solution of nonlinear governing equation.

Concluding Remarks
In this paper, the nonlinear Harry-Dym equation, which describes a couple of nonlinearity and dispersion, has been solved iteratively.It has been shown that Equation (1), which can be transformed to the KdV equation, is an integrable NPDE.Furthermore, it has been shown that the HAM is an efficient method which provides a convenient way to guarantee convergence of the series solutions.The main inferences that can be drawn from this work are summarized as follows: • The HAM will enable us to make an informed choice for auxiliary linear operator and initial approximation which can be employed to obtain the series solution.
• Equation (18) indicates one important role of the auxiliary parameter which was introduced to construct the zeroth-order deformation Equation (4).
• Convergence of the power series (6) largely depends on the range of auxiliary parameter.
• In case of 1 h = − , the HPM is a special case of HAM, as proved by Liao [56].
• By using ( ) , 1 H x t = , influence of the auxiliary function, with the exception of some particular recursive problems (see Ref. [49]) will be neglected.
• Presence of α affects ( ) , u x t considerably when the fractional Harry-Dym equation is utilized.

Figure 1
Figure 1 depicts variations of ( ) , u x t versus x for the Harry-Dym equation calculated by different methods.As shown in this figure, variation of ( ) , u x t

Figure 2 .
Figure 2. Depiction of discrepancy (i.e., error) between the present iterative method and exact solution due to Equation (22).

Figure 3 .
Figure 3. Universal graphs to show the effect of the auxiliary parameter on the convergence of present iterative approach.

Figure 4 .
Figure 4. Universal graphs for predicting the threshold of ( ) , u x t .
took advantage of the HAM to solve linear and nonlinear