Solution of nonlinear higher-index Hessenberg DAEs by Adomian polynomials and differential transform method

The solution of higher-index Hessenberg differential-algebraic equations (DAEs) is of great importance since this type of DAEs often arises in applications. Higher-index DAEs are known to be numerically and analytically difficult to solve. In this paper, we present a new analytical method for the solution of two classes of higher-index Hessenberg DAEs. The method is based on Adomian polynomials and the differential transform method (DTM). First, the DTM is applied to the DAE where the differential transforms of nonlinear terms are calculated using Adomian polynomials. Then, based on the index condition, the resulting recursion system is transformed into a nonsingular linear algebraic system. This system is then solved to obtain the coefficients of the power series solution. The main advantage of the proposed technique is that it does not require an index reduction nor a linearization. Two test problems are solved to demonstrate the effectiveness of the method. In addition, to extend the domain of convergence of the approximate series solution, we propose a post-treatment with Laplace-Padé resummation method.

and incompressible fluids. Unfortunately, these DAEs are known to be difficult to solve, even with numerical methods, due to their complex structure. One reason for this; solutions of higher-index DAEs are constrained for all time by some hidden algebraic equations. As a consequence, initial conditions cannot be prescribed arbitrarily for all solution components as they have to fulfill the constraint equations. Therefore, to start the numerical integration, we need to compute some consistent initial conditions. That is to determine those initial conditions which satisfy all the constraints in the system. Using inconsistent initial conditions or poor estimates can cause the solution of the DAE to drift off the constraints manifold and lead to a non physical solution. Since numerical integration methods have difficulties in solving higher-index DAEs, these problems are usually dealt with by first transforming them to ordinary differential systems (indexzero) or index-one DAEs before applying numerical integration methods. This procedure, known as index-reduction, can be very expensive and may change the properties of the solution of the original problem. Therefore, since important application problems in science and engineering often lead to higher-index DAEs, new techniques are needed to solve these DAEs efficiently.
Over the past decades, significant progress has occurred in the solution of DAEs. Some of these works have focused on the numerical solution and include backward differentiation formula (Brenan 1983), Runge Kutta method (Hairer et al. 1989), pseudospectral method (Hosseini 2005) and finite differences method (Wu and White 2004). One can find other methods for the solution of DAEs like blended implicit methods (Brugnano et al. 2006), implicit Euler (Sand 2002), Chebyshev polynomials (Husein and Jaradat 2008), and arbitrary order Krylov deferred correction methods (Huang et al. 2007).
In recent years, some analytical approximation methods have been developed to solve DAEs. Among such techniques one can find the Adomian decomposition method (ADM) (Hosseini 2006;Celik et al. 2006), the homotopy perturbation method (HPM) (Soltanian et al. 2010;Salehi et al. 2012), the variational iteration method (VIM) (Karta and Celik 2012), the homotopy analysis method (HAM) (Awawdeh et al. 2009), the Padé method (Celik and Bayram 2003) and the differential transform method (DTM) (Benhammouda and Vazquez-Leal 2015; Liu and Song 2007;Ayaz 2004). The ADM, Adomian polynomials and DTM were also applied to solve many other problems. The ADM, for example, was used in computing solutions of algebraic equations (Adomian and Rach 1985;Fatoorehchi et al. , b, 2015Fatoorehchi and Abolghasemi 2014a, b;Fatoorehchi et al. 2015b, d, c). The ADM and Adomian polynomials were applied to various problems in engineering fields (Fatoorehchi et al. 2015f, g, c;Abolghasemi 2015, 2013b). Recently, the DTM was used as a new tool to compute Laplace transforms to solve many problems (Fatoorehchi et al. 2015a;Fatoorehchi and Abolghasemi 2012).
In this work, we present a new procedure for solving nonlinear higher-index Hessenberg DAEs. The method is based on Adomian polynomials (Rach 1984(Rach , 2008Wazwaz 2000;Duan 2010aDuan , b, 2011 and the DTM (Odibat et al. 2010;Lal and Ahlawat 2015;El-Zahar 2013;Fatoorehchi and Abolghasemi 2013a;Gökdoğan et al. 2012;). The DTM is first applied to the DAE where the differential transforms of nonlinear terms are found using Adomian polynomials to obtain a recursion system for the power series coefficients. Based on the index condition, a nonsingular linear recursion system is then derived and solved. It is important to note that the developed procedure does not require an index-reduction nor a linearization. Also it does not depend on complicated tools like perturbation parameters, trial functions, or Lagrangian multipliers as required for perturbation method, HPM or VIM. To enlarge the domain of convergence of the truncated power series, we apply a post-treatment based on Laplace-Padé resummation method Torabi and Yaghoobi 2011;Raftari and Yildirim 2011;Bararnia et al. 2012;George A Baker et al. 1996;Vazquez-Leal et al. 2012;Vazquez-Leal and Guerrero 2014;Khan et al. 2013;.
Two examples of nonlinear higher-index Hessenberg DAEs are solved to demonstrate the effectiveness of the proposed method. Finally, our procedure is straightforward and can be programmed in Maple or Mathematica.
This paper is organized as follows: in "Differential transform method", we review the DTM. Next, in "Padé approximant", "Laplace-Padé resummation method" and "Adomian polynomials and their relation with DTM" we give the basic concepts of Padé approximants, Laplace-Pad é resummation method and Adomian polynomials and their relation with DTM. In "Solution of higher-index Hessenberg DAEs by Adomian polynomials and DTM", we present our analytical method for the solution of nonlinear higher-index Hessenberg DAEs. Then in "Cases study", we apply the developed method to solve two nonlinear higher-index Hessenberg DAEs. Finally, a discussion and a conclusion are given in "Discussion" and "Conclusion", respectively.

Differential transform method
For convenience of the reader, we will review the DTM (Odibat et al. 2010;Lal and Ahlawat 2015;El-Zahar 2013;Fatoorehchi and Abolghasemi 2013a;Gökdoğan et al. 2012) and show how this method is used to solve ordinary differential equations. Definition 2.1 If a function u(t) is analytical with respect to t in the domain of interest, then is the transformed function of u(t).
Definition 2.2 The differential inverse transforms of the set {U k } n k=0 is defined by From the above definitions, it is easy to see that the concept of the DTM is obtained from the power series expansion. To illustrate the application of the DTM to solve ordinary differential equations, we consider the nonlinear equation where f (u(t), t) is a nonlinear smooth function. Equation (4) is supplied with some initial condition DTM establishes that the solution of (4) can be written as where U 0 , U 1 , U 2 , . . . are unknowns to be determined by DTM.
Applying the DTM to the initial condition (5) and equation ( 4) respectively, we obtain the transformed initial condition and the recursion equation where F U 0 , . . . , U k−1 , k − 1 is the differential transforms of f (u(t), t).
Using (7) and (8), we determine the unknowns U k , k = 0, 1, 2, . . . Then, the differential inverse transformation of the set of values {U k } m k=0 gives the approximate solution where m is the approximation order of the solution. The exact solution of problem (4-5) is then given by (6). If U k and V k are the differential transforms of u(t) and v(t) respectively, then the main operations of DTM are shown in Table 1.
The process of the DTM can be described as: 1. Apply the differential transform to initial condition (5). 2. Apply the differential transform to the differential equation ( 4) to obtain a recursion equation for the unknowns U 0 , U 1 , U 2 , . . . 3. Use the transformed initial condition (7) and the recursion equation (8) to determine the unknowns U 0 , U 1 , U 2 , . . . 4. Use the differential inverse transform formula (9) to obtain an approximate solution for initial-value problem (4-5). (4) The solutions series obtained from DTM may have limited regions of convergence. Therefore, we propose to apply the Laplace-Padé resummation method to DTM truncated series to enlarge the convergence region as depicted in the next sections.

Given an analytical function u(t) with Maclaurin's expansion
The Padé approximant to u(t) of order [L, M] which we denote by [L/M] u (t) is defined by George A Baker et al. (1996) where we considered q 0 = 1, and the numerator and denominator have no common factors.
The numerator and the denominator in (11) are constructed so that u(t) and [L/M] u (t) and their derivatives agree at t = 0 up to L + M. That is From (12), we have From (13), we get the following algebraic linear systems From (14), we calculate first all the coefficients q n , 1 ≤ n ≤ M. Then, we determine the coefficients p n , 0 ≤ n ≤ L from (15). Note that for a fixed value of L + M + 1, the error (12) is smallest when the numerator and denominator of (11) have the same degree or when the numerator has degree one higher than the denominator.

Laplace-Padé resummation method
Several approximate methods provide power series solutions (polynomial). Nevertheless, sometimes, this type of solutions lack large domains of convergence. Therefore, Laplace-Padé resummation method is used in literature to enlarge the domain of convergence of solutions or to find the exact solutions.
The Laplace-Padé method can be summarized as follows: 1. First, Laplace transformation is applied to power series (9). 2. Next, s is substituted by 1/t in the resulting equation.
3. After that, we convert the transformed series into a meromorphic function by forming its Padé approximant of order [N/M]. N and M are arbitrarily chosen, but they should be smaller than the order of the power series. In this step, the Padé approximant extends the domain of the truncated series solution to obtain better accuracy and convergence. 4. Then, t is substituted by 1/s. 5. Finally, by using the inverse Laplace s transformation, we obtain the exact or an approximate solution.

Adomian polynomials and their relation with DTM
In this section, we briefly review the Adomian polynomials and their relation with the DTM. Usually a nonlinear term N(u) in a differential equation is decomposed in terms of Adomian polynomials A n (Rach 2008(Rach , 1984Wazwaz 2000;Duan 2010aDuan , b, 2011 as where A n are generated for all forms of nonlinearity from and where u n (t), n = 0, 1, 2, . . . denote the components used in the expansion There are several algorithms to compute Adomian polynomials but recently a convenient recursion to calculate Adomian polynomials for the m-variable case is proposed in (Duan 2011) Also an extension of the differential transform to nonlinear terms of any type, known as the improved DTM, was given in (Fatoorehchi and Abolghasemi 2013a, 2014b) using Adomian polynomials In the coming sections, we make use of (19) and (20 ) to show how to solve nonlinear higher-index Hessenberg DAEs.

Solution of higher-index Hessenberg DAEs by Adomian polynomials and DTM
In this section, we present our method for solving nonlinear higher-index Hessenberg differential-algebraic equations (DAEs). The technique is based on Adomian polynomials and the differential transform method (DTM). To solve the DAE, we first apply the DTM to it, where Adomian polynomials are used to compute the differential transforms of the nonlinear terms. The resulting recursion equations are rearranged in a nonsingular linear algebraic system for the coefficients of the power series solution. Two classes of nonlinear higher-index Hessenberg DAEs are solved.

Higher-index nonlinear Hessenberg DAEs
The first class of higher-index Hessenberg DAEs we consider here is The DAE is supplied with some consistent initial conditions η i are given constants. System (21-22) has index (m + 1) if the product of the Jacobians is nonsingular for t ≥ 0.

(19)
An important subclass of system (21-22) consists of those DAEs arising from the simulation of constrained mechanical multibody systems. Such DAEs have the form where u(t) is the vector of generalized coordinates, ü(t) is the vector that contains the system accelerations, ∂g/∂u is the Jacobian of g, v(t) is the Lagrange multipliers vector and f (u(t)) is the generalized forces vector.
A standard assumption for these DAEs is the full rank condition which means that the constraint equations are linearly independent. If condition (27) is satisfied then , then using (19), the Adomian polynomials F j k , j = 1, . . . , n u , k = 0, 1, 2, . . . for the (n u + n v )-variable function f j (u, v) are given by where U i,l and V i,l are the differential transforms of u i and v i .
Equation (30) can be written as In vector form, we have In a similar manner, let g(u) = g 1 (u), g 2 (u), . . . , g n v (u) T then the Adomian polynomials G j k , j = 1, . . . , n v , k = 0, 1, 2, . . . for the n u -variable function g j (u) are given by In vector form, we have To solve DAE (21-22), we apply the DTM to get and where U k is the differential transform of u(t) and α = k(k − 1) . . . (k + 1 − m).

System (40) can be decomposed as
Since condition (24) holds, then the first equation of (43) can be solved uniquely for V k−m . Then using the second equation of (43), we can determine U k . Therefore, an approximate analytical solution is given by

Index-three nonlinear Hessenberg DAEs
The second class of higher-index nonlinear Hessenberg DAEs we consider here is The DAE is supplied with some consistent initial conditions System (45) is index-three if the product of the Jacobians is nonsingular for t ≥ 0.
Let us assume that f, g and h are sufficiently smooth and that the Jacobian ∂g/∂u has full row rank [i.e. rank ∂g/∂u = n w ] for t ≥ 0. Let T then the Adomian polynomials F j k , j = 1, . . . , n u , k = 0, 1, 2, . . . for the (n u + n v )-variable function f j (u, v) are given by Equation (49) can be written as

In vector form, we have
where T then the Adomian polynomials H j k , j = 1, . . . , n v , k = 0, 1, 2, . . . for the (n u + n v + n w )-variable function h j (u, v, w) are given by Equation (54)  . In a similar manner, let g(u) = g 1 (u), g 2 (u), . . . , g n v (u) T then the Adomian polynomials G j k , j = 1, . . . , n v , k = 0, 1, 2, . . . for the n u -variable function g j (u) are given by To solve DAE (45-46), we apply the DTM to get and where U k , V k and W k are the differential transforms of u(t), v(t) and w(t) respectively.
From the (61), we finally come to the linear recursion system where System (63) can be decomposed as Since condition (47) holds, then the first equation of (65) can solved uniquely for W k−2 . Then V k−1 is obtained from the second equation of (65). Last, the unknown U k is obtained from the third equation of (65). Then, an approximate analytical solution is given by

Cases study
In this section, we will demonstrate the effectiveness of proposed technique through two nonlinear higher-index Hessenberg DAEs.

Example 1
Consider the following nonlinear index-three Hessenberg DAE describing the constrained motion of a particle to a circular track System (67) is supplied with the following (consistent) initial conditions Note that no initial condition v(0) is given to the variable v(t) as v(0) is pre-determined by the DAE and initial conditions (68). System (67) is index-three since three time differentiations of the algebraic equation (third equation) of (67) will lead to an ordinary differential equation for v(t). As a consequence, this DAE system is difficult to solve numerically due to numerical instabilities. Therefore, to solve (67-68), we apply the DTM to (67) and get the recursion where the differential transform of the nonlinear terms u 3 i (t), i = 1, 2 are replaced by the Adomian polynomials Then applying the DTM to initial conditions (68), we get For k = 0 and k = 1, the third equation of (69) gives which are satisfied by the transformed initial conditions (70).

Example 2
Consider the following nonlinear index-three Hessenberg DAE where System (79) is supplied with the following (consistent) initial conditions Note that no initial condition w(0) is given to the variable w(t) as w(0) is pre-determined by the DAE and initial conditions (80). System (79) is index-three since three time differentiations of the algebraic equation (fifth equation) of (79) will lead to an ordinary differential equation for w(t). As a consequence, this DAE system is difficult to solve numerically due to numerical instabilities.
To solve (79-80), we first expand ϕ 1 (t) and ϕ 2 (t) in Taylor series Then, we apply the DTM to (79) and get the recursion where i,k is the differential transform of ϕ i (t), for i = 1, 2, 3 and where the differential transform of the nonlinear terms e u i , i = 1, 2 are replaced by the Adomian polynomials A i k (78) u 1 (t) = cos t, u 2 (t) = sin t, v(t) = 1 + sin 2t, kU 1,k = 2V 1,k−1 , kU 2,k = 2V 2,k−1 , Then, we apply the DTM to initial conditions (80), to get Using the first two equations of (82) with k = 1 and (83), we get For k = 0 and k = 1, the last equation of (82) gives which are satisfied by (83) and (84).
Therefore, system (82) reduces to the following nonsingular linear algebraic system for the unknowns U 1,k , U 2,k , V 1,k−1 , V 2,k−1 and W k−2 Adding the third and the fourth equations and using the last equation, we obtain W k−2 . Now replacing W k−2 by its expression in third and fourth equations, we get U 1,k and U 2,k . Last, we use the first and second equations to obtain V 1,k−1 and V 2,k−1 . Following this procedure and using (83) and (84), we obtain the approximations which are the first terms of the Taylor series of the exact solutions of DAE initial-value problem (79-80).

Discussion
Higher-index differential-algebraic equations (DAEs) still require new numerical and analytical methods to solve them efficiently. Such problems are known to be difficult to solve both numerically and analytically. In this paper, we introduced a new analytical method to solve nonlinear higher-index Hessenberg DAEs. The method is based on Adomian polynomials and the differential transform method (DTM). Two classes of nonlinear higher-index Hessenberg DAEs were treated by this method. The method has successfully handled these two classes of DAEs without the need for a preprocessing step of index-reduction. The method transformed the DAEs into easily solvable linear algebraic systems for the coefficient of the power series solution. For each class, one test problem was solved. The examples show that Adomian polynomials combined with the DTM are powerful tools to obtain the exact solutions or approximate solutions of nonlinear higher-index Hessenberg DAEs. To improve the power series solution, a Laplace-Padé post-treatement is applied to the truncated series leading to the exact solution.

Conclusion
This work presents the analytical solution of two classes of nonlinear higher-index Hessenberg DAEs using Adomian polynomials and the DTM. Procedures for solving these two classes of DAEs are presented. For each class, the technique was tested on one nonlinear higher-index Hessenberg problem. The results obtained show that the method can be applied to solve nonlinear higher-index Hessenberg DAEs efficiently obtaining the exact solution or an approximate solution. On the one hand, it is important to note that these types of DAEs are difficult to solve both numerically and analytically. On the other hand, the presented technique based on Adomian polynomials and the DTM in combination with Laplace-Padé resummation method was able to obtain the exact solution of nonlinear higher-index Hessenberg DAEs. The use of Adomian polynomials allowed us to obtain an algorithm for the method and also to compute the differential transforms of highly nonlinear terms. The technique is based on a straightforward procedure that can be programmed in Maple or Mathematica to simulate large problems. Finally, future work is needed to apply the proposed technique to higher-index partial differentialalgebraic equations and other nonlinear higher-index DAEs. Our method can be combined with the multi-stage DTM to calculate accurate approximate solutions to these problems.

Competing interests
The author declares that there is no conflict of interests regarding the publication of this paper.