# Preface

This section discusses a simplified version of the Adomian decomposition method first concept of which was proposed by Randolf Rach in 1989 that was crystallized later in a paper published with his colleagues G. Adomian and R.E. Meyers. That is way this technique is frequently referred to as the Rach--Adomian--Meyers modified decomposition method (MDM for short). Initially, this method was applied to power series expansions, which was based on the nonlinear transformation of series by the Adomian--Rach Theorem. Similar to the Runge--Kutta methods or spline method, the MDM can be implemented in numerical integration of differential equations by one-step methods. In case of polynomials or power series, it shows the advantage in speed and accuracy of calculations when at each step the Adomian decomposition method allows one to perform explicit evaluations.

# Modified Decomposition Method

The modified decomposition method (MDM) is both a decelerated Adomian decomposition method and a developed power series method that nicely and accurately treats any analytic nonlinearity in differential equations. In case of autonomous homogeneous or nonhomogeneous equations, the method is seems to be superior over other methods. Therefore, we start demonstration of this method in initial value problems for first order differential equations of the form

$L \left[ y(t) \right] = N\left[ y(t) \right] + g(t) , \qquad y(t_0 ) = y_0 ,$
where $$L \left[ y(t) \right] = y' (t) + a (t)\, y (t) = {\text d}y/{\text d}t + a (t)\, y (t)$$ is the linear differential operator, $$R \left[ y(t) \right] = b (t)\,y(t) ,$$ is a linear operator, and $$N \left[ y(t) \right]$$ is a nonlinear term, g is the system input and y is the system output. The operator can be rewritten as $$L = \texttt{D} + a(t)\,\texttt{T}$$ , where $$\texttt{D} = {\text d}/{\text d}t$$ is the derivative operator and $$\texttt{I}$$ is the identical operator. Since the unbounded derivative operator has a one-dimensional kernel (the set of functions that are annihilated by the derivative operator), its inverse, calling antiderivative, is a multivariable operator depending on some constant. Similar property is valid for L. In order to knock out a single inverse operator L-1, we consider consider L on a space of functions that satisfy an initial condition; more over, it is assumed that an explicit formula is available. Furthermore, we emphasize that the choice for L and concomitantly its inverse L-1 are determined by the particular equation to be solved. Hence, the choice for L is nonunique, e.g. for cases of differential equations with singular coefficients, a different form for the linear operator L can be chosen. As a rule, the highest derivative is taken as the linear operator L. Later we consider more general equation including coefficients depending on independent variable.

The MDM and its modifications and generalizations have been extensively applied in physics, chemistry, mechanics, hydrology, engineering, economics, biology, epidemiology, etc. Some references can be found in the following link.

# Complete List of Publications using the Modified Decomposition Method

The Rach--Adomian--Meyers MDM decomposes the solution into a power series

$y(t) = \sum_{n\ge 0} c_n \left( t- t_0 \right)^n$
and decomposes the nonlinearity into a power series according to the Adomian--Rach theorem
$N \left[ y(t) \right] = N \left[ \sum_{n\ge 0} c_n \left( t- t_0 \right)^n \right] = \sum_{n\ge 0} A_n \left( c_0 , c_1 , \ldots , c_n \right) \left( t- t_0 \right)^n ,$
where An are the Adomian polynomials in terms of solution coefficients. Of course, the above formula is nothing more than Taylor's series expansion of the nonlinear function N[y]. However, the main difference of the MDM expansion from the power series is that in the former the coefficients depend recursively on the previous coefficients whereas Taylor's coefficients all depend on the infinitesimal behavior of the slope function, but not on previously found coefficients. It is convenient to make a shift in independent variable x = t - t0 to transfer the initial point from t0 into the origin:
$y(x) = y_0 + \sum_{n\ge 1} c_n \, x^n \qquad\mbox{and} \qquad N \left[ y(x) \right] = \sum_{n\ge 0} A_n \left( c_0 , c_1 , \ldots , c_n \right) x^n .$
The Adomian polynomials are evaluated in a similar way as in the ADM
$A_n = \frac{1}{n!} \, \frac{\partial^n}{\partial \lambda^n} \left( N\left[ \sum_{n\ge 0} u_n \lambda^n \right] \right)_{\lambda =0} , \qquad n=0,1,2,\ldots ;$
but with substitution cn instead of un:
$A_n \left( u_0 , u_1 , \ldots , u_n \right) = x^n A_n \left( c_0 , c_1 , \ldots , c_n \right) , \qquad n=0,1,2,\ldots .$
Actually, Mathematica can find any Adomian polynomail for polynomial nonlinearities with a dedicated command Coefficient or CoefficientList; so sometimes the above formula is not needed when there is an access to the computer algebra system. Recall that if N is not a polynomail, then Mathematica cannot provide you an answer, and you have to use a special script to determine coefficients.

CoefficientList[(c0 + c1 x + c2 *x^2 + c3 x^3 + c4 x^4 + c5 x^5)^2 , x]
{c0^2, 2 c0 c1, c1^2 + 2 c0 c2, 2 c1 c2 + 2 c0 c3, c2^2 + 2 c1 c3 + 2 c0 c4, 2 c2 c3 + 2 c1 c4 + 2 c0 c5, c3^2 + 2 c2 c4 + 2 c1 c5, 2 c3 c4 + 2 c2 c5, c4^2 + 2 c3 c5, 2 c4 c5, c5^2}

For example, if N[y] = y², we have $$A_n = \sum_{k=0}^n u_k u_{n-k} .$$ Therefore,

A0(u0) = u0² A0(c0) = c0²
A1(u0, u1) =2 u0 u1 A1(c0, c1) =2 c0 c1
A2(u0, u1, u2) = u1² + 2 u0 u2 A2(c0, c1, c2) = c1² + 2 c0 c2
A3(u0, u1, u2, u3) = 2 u1u2 + 2 u0 u3 A3(c0, c1, c2, c3) = 2 c1c2 + 2 c0 c3
Coefficient[(u0 + s*u1 + s^2 *u2 + u3*s^3 + s^4 *u4)^2, s, 3]
2 u1 u2 + 2 u0 u3
Mathematica code for evaluating Adomian polynomials (Jun-Sheng Duan, An efficient algorithm for the multivariable Adomian polynomials, Applied Mathematics and Computation, 2010, 217, pp. 2456--2467):
Subscript[A, 0] = f[Subscript[u, 0]];
Z = Table[0, {i, 1, n}, {j, 1, i}];
Do[Z[[m, 1]] = Subscript[u, m], {m, 1, n}];
For[m = 2, m <= n, m++,
For[k = 2, k <= m, k++,
Z[[m, k]] = Expand[Subscript[u, 1]*Z[[m - 1, k - 1]]];
If[Head[Z[[m, k]]] === Plus, Z[[m, k]] = Map[#/Exponent[#, Subscript[u, 1]]&, Z[[m, k]]], Z[[m, k]] = Map[#/Exponent[#, Subscript[u, 1]]&, Z[[m, k]],{0}]]];
For[k = 2, k <= Floor[m/2], k++, Z[[m, k]] = Z[[m, k]] + (Z[[m - k, k]]/. Subscript[u, s_] -> Subscript[u, s + 1])]];
dir = Table[D[f[Subscript[u, 0]], {Subscript[u, 0], k}], {k, 1, n}];
Do[Subscript[A, m] = Take[dir, m].Z[[m]], {m, 1, n}];
Table[Subscript[A, m], {m, 0, n}]]
Here is its another version (Jun-Sheng Duan):
pol = Series[f[Sum[Subscript[u, i]*n[Lambda] ^ i,
{i, 0, n}]], {n[Lambda], 0, n}];
For[i = 0, i<=n, i++, Subscript[A,i]=Collect[ExpandAll[ SeriesCoefficient[pol,i]],Derivative[_][f][_]]];
Table[Subscript[A, i], {i, 0, n}]]
Before we go to more sofisticated examples, we demonstrate the MDM on the following familiar one. Let $$u = e^x = \sum_{n\ge 0} \frac{x^n}{n!} = \sum_{n\ge 0} c_n x^n .$$ Raising it to the third power, we get
$\left( e^x \right)^3 = \sum_{n\ge 0} A_n x^n = e^{3x} = \sum_{k\ge 0} \frac{(3x)^k}{k!} = \sum_{k\ge 0} \left( \frac{3^k}{k!} \right) x^k \qquad \Longrightarrow \qquad A_n = \frac{3^n}{n!} .$
Now we evaluate the Adomian polynomials using  CoefficientList
CoefficientList[(1 + 1 x + 1/2 *x^2 + 1/3! x^3 + 1/4! x^4 + 1/5! x^5)^3 , x]
{1, 3, 9/2, 9/2, 27/8, 81/40, 121/120, 17/40, 49/320, 409/8640, \ 361/28800, 9/3200, 181/345600, 1/12800, 1/115200, 1/1728000}
CoefficientList[(c0 + c1 x + c2 *x^2 + c3 x^3 + c4 x^4 + c5 x^5)^3 , x]
\begin{align*} A_0 \left( u_0 \right) &= A_0 \left( c_0 \right) = c_0^3 = 1, \\ A_1 \left( u_0 , u_1 \right) &= 3\,c_0^2 c_1 = 3 , \\ A_2 \left( u_0 , u_1 , u_2 \right) &= 3\,c_0 c_1^2 + 3\,c_0^2 c_2 = 9/2 , \\ A_3 \left( u_0 , u_1 , u_2 , u_3 \right) &= c_1^3 + 6\. c_0 c_1 c_2 + 3\, c_0^2 c_3 = 9/2 , \\ A_4 \left( u_0 , u_1 , u_2 , u_3 , u_4 \right) &= 3\, c_1^2 c_2 + 3\, c_0 c_2^2 + 6\, c_0 c_1 c_3 + 3\, c_0^2 c_4 = 27/8 , \end{align*}
and so on.

Example: Let us consider the arccotangent function on the interval |x| < 1

$\mbox{arccot} (x) = \frac{\pi}{2} - \arctan (x) = \frac{\pi}{2} -x + \frac{x^3}{3} - \frac{x^5}{5} + \frac{x^7}{7} - \cdots = \frac{\pi}{2} + \sum_{k\ge 1} (-1)^k \frac{x^{2k+1}}{2k+1} , \qquad |x| < 1.$
This means that we have a Maclaurin expansion of this function that converges on the in the open interval |x| < 1:
$\mbox{arccot} (x) = \sum_{k\ge 0} c_k x^k \qquad\mbox{with} \quad c_0 = \frac{\pi}{2} , \quad c_{2k} =0 \quad\mbox{for} \quad k=1,2,\ldots , \quad c_{2k+1} = \frac{(-1)^k}{2k+1}\quad\mbox{for} \quad k=1,2,\ldots .$
Assuming[x > 0, Series[ArcCot[x], {x, 0, 8}]]
$$\frac{\pi}{2} -x + \frac{x^3}{3} - \frac{x^5}{5} + \frac{x^7}{7} + O\left[ x \right]^9$$
Suppose we want to find the Maclaurin seris for 1/(arccot(x))³; so we first use Mathematica to obtain its correct expansion:
f[y_] = 1/(ArcCot[y])^3;
Assuming[x > 0, Series[f[x], {x, 0, 9}]]
which yields
$\frac{1}{\left( \mbox{arccot} (x) \right)^3} = \frac{8}{\pi^3} + \frac{48}{\pi^4}\,x + \frac{192}{\pi^5}\, x^2 - \frac{16 \left( \pi^2 - 40 \right)}{\pi^6}\, x^3 - \frac{128 \left( \pi^2 -15 \right)}{\pi^7}\, x^4 + \frac{16 \left( 9\pi^6 -784 \pi^4 + 11760 \pi^2 -48384 \right)}{5\, \pi^8}\, x^5 + \cdots$
Choosing f(y) = 1/y³, we expand the composition of two functions f and arccot into power series using Adomian's polynomials:
$\frac{1}{\left( \mbox{arccot} (x) \right)^3} = \left( \frac{2}{\pi} \right)^3 + \sum_{k\ge 1} A_k \left( c_0 , c_1 , \ldots , c_k \right) x^k ,$
where the first term is obviously
$A_0 (c_0 ) = f(c_0 ) = \frac{1}{c_0^3} , \qquad \mbox{with} \quad c_0 = \frac{\pi}{2} .$
Using the derivative of the arccotangent function, $${\text d}\,\mbox{arccot}(x)/{\text d}x = -\left( 1 + x^2 \right)^{-1} ,$$ we calculate the next Adomian's polynomial:
$A_1 = u_1 f'(u_0 ) = -\frac{3\,u_1}{u_0^4} = x\,c_1 f' (c_0 ) = -\frac{3\,c_1}{c_0^4} = \frac{48}{\pi^4} ,$
where c1 = -1, which follows from arccotangent Maclaurin power series expansion. The second Adomian polynomail is
$A_2 = u_2 f' (u_0 ) + \frac{1}{2}\, u_1^2 f'' (u_0 ) = - \frac{3\,u_2}{u_0^4} + \frac{6\, u_1^2}{u_0^5}$
Therefore,
$A_2 \left( c_0 , c_1 , c_2 \right) = \frac{6\, c_1^2}{c_0^5} - \frac{3\,c_2}{c_0^4} = \frac{6}{c_0^5} = \frac{6*2^5}{\pi^5} = \frac{192}{\pi^5}$
because c2 = 0 and c1 = -1. Next term is
$A_3 = u_3 f'(u_0 ) + \frac{1}{6}\,u_1^3 f''' (u_0 ) + u_1 u_2 f'' (u_0 ) = -u_3 \frac{3}{u_0^4} - \frac{60}{6}\,u_1^3 \frac{1}{u_0^6} + u_1 u_2 \,\frac{12}{u_0^5} .$
Using values c0 = π/2, c1 = -1, c2 = 0, and c3 = 1/3, we obtain
$A_2 \left( c_0 , c_1 , c_2 , c_3 \right) = A_2 \left( \frac{\pi}{2} , -1 , 0 , \frac{1}{3} \right) = -c_3 \frac{3}{c_0^4} + \frac{10}{c_0^6} = \frac{10}{c_0^6} - \frac{1}{c_0^4} = \frac{10-c_0^2}{c_0^6} = 2^4 \frac{40- \pi^2}{\pi^6} .$
■

Once we expand into Maclauirin power series every term in the given differential equation (assuming that shifting of independent variable was done)

$\frac{{\text d}y}{{\text d}x} + a(x)\,y(x) = b(x)\,y(x) + N \left[ y(x) \right] + g(x)$
we obtain
$\sum_{k\ge 1} c_k \,k \,x^{k-1} + \left( \sum_{k\ge 0} a_k x^k \right) \left( \sum_{k\ge 0} c_k x^k \right) = \left( \sum_{k\ge 0} b_k x^k \right) \left( \sum_{k\ge 0} c_k x^k \right) + \sum_{k\ge 0} A_k \left( c_0 , \ldots , c_k \right) x^k + \sum_{k\ge 0} g_k x^k ,$
where coefficients ck are unknown yet. Using the convolution rule for power series (also known as Cauchy product), we get
$\sum_{k\ge 1} c_k \,k \,x^{k-1} + \sum_{k\ge 0} \left( \sum_{i=0}^k a_{k-i} c_i \right) x^k = \sum_{k\ge 0} \left( \sum_{i=0}^k b_{k-i} c_i \right) x^k + \sum_{k\ge 0} A_k \left( c_0 , \ldots , c_k \right) x^k + \sum_{k\ge 0} g_k x^k .$
Then equating coefficients of like powers of x, we derive a recurrence relation for determination of unknown coefficients ck.

$c_{k+1} \left( k+1 \right) + \sum_{i=0}^k a_{k-i} c_i = \sum_{i=0}^k b_{k-i} c_i + A_k \left( c_0 , \ldots , c_k \right) + g_k , \quad c_0 = y(0) , \qquad k=0,1,2,\ldots .$

Since we are working with holomorphic functions, it is natural to assume that the input function g is represented by a convergent power series

$g(x) = \sum_{n\ge 0} g_n \frac{x^n}{n!} .$
Substituting all these series into the differential equation, we obtain
\begin{align*} L \left[ y(x) \right] &= \texttt{D} \sum_{n\ge 1} c_n x^n + a(x)\,\sum_{n\ge 0} c_n x^n = \sum_{n\ge 1} n\,c_n x^{n-1} + \left( \sum_{n\ge 0} \alpha_n x^n \right) \left( \sum_{n\ge 0} c_n x^n \right) , \\ R \left[ y(x) \right] &= b(x)\,y(x) = \left( \sum_{n\ge 0} \beta_n x^n \right) \left( \sum_{n\ge 0} c_n x^n \right) , \\ N \left[ y(x) \right] &= \sum_{n\ge 0} A_n x^n . \end{align*}
Since the multiplication of two series leads to their convolution,
$\left( \sum_{n\ge 0} \a_n x^n \right) \left( \sum_{n\ge 0} b_n x^n \right) = \sum_{n\ge 0} \left( \sum_{k=0}^n a_k b_{n-k} \right) x^n$
we get the equation for power series:
$\sum_{n\ge 1} n\,c_n x^{n-1} + \sum_{n\ge 0} \left( \sum_{k=0}^n \alpha_k c_{n-k} \right) x^n = \sum_{n\ge 0} \left( \sum_{k=0}^n \beta_k c_{n-k} \right) x^n + \sum_{n\ge 0} A_n x^n + \sum_{n\ge 0} g_n \frac{x^n}{n!} ,$
with
$c_0 = y_0 .$

Example: Consider the initial value problem

$\frac{{\text d}y}{{\text d}x} + 2x\, y(x) = 2 + 2\,x^2 + \frac{y^2}{2} , \qquad y(0) =1.$
The given problem has the exact solution
$y(x) = 2x + \frac{1}{1-x/2} = 1 + \frac{5x}{2} + \frac{x^2}{4} + \frac{x^3}{8} + \left( \frac{x}{2} \right)^4 + \left( \frac{x}{2} \right)^5 + \left( \frac{x}{2} \right)^6 + \cdots + \left( \frac{x}{2} \right)^n + \cdots .$
which we verify with Mathematica:
y[x_] = 2*x + 1/(1 - x/2);
Simplify[D[y[x], x] + 2*x*y[x] - 2 - 2*x^2 - (1/2)*(y[x])^2]
0
Series[y[x], {x, 0, 10}]
Let L be the differential operator
$L\left[ x, \texttt{D} \right] = e^{-x^2} \texttt{D} e^{x^2} = \texttt{D} + 2x\,\texttt{I} \qquad \Longrightarrow \qquad L^{-1}\left[ x, \texttt{D} \right] f(x) = e^{-x^2} + e^{-x^2} \int_0^x e^{s^2} f(s)\,{\text d}s ,$
and $$N\left[ y \right] = y^2$$ be the nonlinear operator. Then we rewrite the given problem in the operator form:
$L\left[ x, \texttt{D} \right] y = 2 + 2\,x^2 + \frac{1}{2}\,N\left[ y \right]$
Assuming that the solution can be decomposed into convergent power series
$y(x) = 1 + \sum_{n\ge 1} c_n x^n ,$
we expand the nonlinear term into the series :
$N \left[ y \right] = y^2 = N \left[ 1 + \sum_{n\ge 1} c_n x^n \right] = \sum_{n\ge 0} A_n \left( 1, c_1 , \ldots , c_n \right) x^n ,$
where An are Adomian's polynomials. Substituting these series into the given differential equation, we get
\begin{align*} L\left[ x, \texttt{D} \right] y(x) &= \sum_{n\ge 1} n\,c_n x^{n-1} + 2x + 2 \sum_{n\ge 1} c_n x^{n+1} = 2 + 2\,x^2 + \frac{1}{2} \sum_{n\ge 0} A_n x^n . \end{align*}
The first Adomian polynomails should be
$A_0 = 1 = N\left[ 1 \right]$
because of the initial condition y(0) = 1 = c0. Now we compare coefficients of like powers and obtain the sequence of recurrence relations:
\begin{align*} c_1 &= 2 + \frac{1}{2}\, A_0 = \frac{5}{2} , \\ 2\,c_2 +2 &= \frac{1}{2}\, A_1 = c_1 c_0 \qquad \Longrightarrow \qquad c_2 = \frac{1}{4} , \\ 3\, c_3 + 2\, c_1 &= 2 + \frac{1}{2}\, A_2 = 2 + c_0 c_2 + \frac{1}{2}\, c_2^2 \qquad \Longrightarrow \qquad c_3 = \frac{1}{8} , \\ \left( n+1 \right) c_{n+1} + 2\,c_{n-1} &= \frac{1}{2}\, A_n , \qquad n= 3,4, \ldots . \end{align*}
Therefore, we get a recurrence:
$c_{n+1} = \frac{1}{2(n+1)} \, A_n \left( 1, c_1, c_2 , \ldots , c_n \right) - \frac{2}{n+1}\, c_{n-1} = \frac{1}{2^{n+1}} , \qquad n=3,4,\ldots .$
This formulas can be verified by induction.    ■

Example: Consider the first-order nonlinear differential equation that contains arbitrary holomorphic function

$\frac{{\text d} y}{{\text d}x} +2x\,f \left( y(x) \right) = 3x, \qquad y(0) =1,$
where f is a known function having required number of derivatives (generally speaking, all derivatives). Assuming that the given differential equation has a convergent Maclaurin series solution
$y(x) = \sum_{n\ge 0} a_n x^n ,$
we expand the nonlinear term into Adomian's series
$f \left( y(x) \right) = \sum_{n\ge 0} A_n \left( a_0 , a_1 , \ldots , a_n \right) x^n = \sum_{n\ge 0} A_n x^n .$
Then the derivative becomes
$y' (x) = \sum_{n\ge 0} (n+1)\, a_{n+1} x^n .$
Similar, we expand the complete nonlinear term:
$x\,f \left( y(x) \right) = \sum_{n\ge 0} A_n \left( a_0 , a_1 , \ldots , a_n \right) x^{n+1} = A_0 x + \sum_{n\ge 0} A_{n+1} x^{n+2} ,$
where An are Adomian's polynomials. Substituting these identities into the given differential equation, we obtain
$a_1 + 2\,a_2 x + \sum_{n\ge 0} \left( n+3 \right) a_{n+3} x^{n+2} = 3x - 2\,A_0 x -2\,\sum_{n\ge 0} A_{n+1} x^{n+2} .$
From the initial condition, we find immediately that a0 = 1, and for all other coefficients, we get the recurrence relation by equating the coefficients of like powers of x:
\begin{align*} a_1 &= 0 , \\ 2\,a_2 &= 3 - 2\, A_0 = 3 - 2\, f(a_0 ) \qquad \Longrightarrow \qquad a_2 = \frac{1}{2} \left( 3 - 2\,f(1) \right) , \\ a_{n+3} &= - \frac{2}{n+3} \, A_{n+1} , \qquad n= 0,1,2,\ldots . \end{align*}
The truncated sum
$\phi_m (x) = \sum_{k=0}^m a_k x^k$
gives a polynomial approximation to the solution.    ■

Example: Consider the following Riccati differential equation:

$\frac{{\text d} y}{{\text d}x} = x^2 + y^2 (x), \qquad y(x_0 ) = y_0 .$
We rewrite the given differential equation in operator form:
$L\left[ y \right] = f(y) + g(x) ,$
where
$L\left[ y \right] = \texttt{D}\,y(x) = \frac{{\text d} y}{{\text d}x} , \qquad g(x) = x^2 , \quad f(y) = y^2 (x) .$
The inverse of the unbounded differential operator L is
$L^{-1} \left[ u \right] = u(t_0 ) + \int_{x_0}^x u(t)\,{\text d}t .$
Application of the bounded operator to the given differential equation yields the functional equation for which the fixed point theorem can be applied:
$y(x) = L^{-1} \left[ g(x) \right] + L^{-1} \left[ f(y) \right] .$
The Adomian decomposition method (ADM) states that the dependent variable y(x) and the nonlinear term f(y) = y² should be written as the following infinite series
\begin{align*} y(x) &= \sum_{n\ge 0} u_k (x) , \\ f(y(x)) &= \sum_{n\ge 0} A_k \left( u_1 , u_2 , \ldots , u_k \right) , \end{align*}
where Ak is the k-th Adomian polynomial.    ■

Example: Consider the initial value problem

$y' = e^y + \tan (2x) - \cos (2x), \qquad y(0) =0;$
that has the true solution $$y(x) = \ln \left( \cos 2x \right) .$$ With Mathematica, we find its power series expansion is
Series[Log[Cos[2*x]], {x, 0, 12}]
$y(x) = -2x - \frac{4}{3}\, x^4 - frac{64}{45}\, x^6 - \frac{544}{315}\, x^8 - \frac{1415168}{46777}\, x^{10} - \frac{1415168}{467775}\, x^{12} - \cdots .$
According to the standard ADM, we assume that solution is represented by infinite series
$y(x) = u_0 (x) + \sum_{k\ge 1} u_k (x) ,$
where
$u_0 (x) = \tan (2x) - \cos (2x) ,$
and all other terms are determined recursively by solving the initial value problems:
$u'_{m+1} (x) = A_m \left( u_0 , u_1 , \ldots u_m \right) , \qquad u_{m+1} (0) = 0, \qquad m=0,1,2,\ldots .$
Since the nonlinear operator in the given equation is $$N[y] = e^y ,$$ we get the first Adomian polynomial directly
$A_0 (u_0 ) = N[ u_0 ] = e^{u_0} = e{\tan (2x) - \cos (2x)} .$
■

General first order differential equation

We consider the following initial value problem for the nonlinear differential equation with holomorphic input function g(x) and holomorphic coefficients
$\frac{{\text d}y}{{\text d}x} = a(x)\,y(x) + b(x)\,N \left[ y(x) \right] + g(x) , \qquad y(0) = y_0 .$
We assume that this problem has a Maclaurin power series solution
$y (x) = c_0 (x) + \sum_{k\ge 1} c_k \,x^k , \qquad \mbox{with} \quad c_0 = y(0) = y_0 .$
Since all known functions in the above differential equation are holomorphic, we have
\begin{align*} a(x) &= \sum_{k\ge 0} a_k x^k , \\ b(x) &= \sum_{k\ge 0} b_k x^k , \\ g(x) &= \sum_{k\ge 0} g_k x^k , \end{align*}

Example: Consider the first-order nonlinear differential equation with a sinusoidal nonlinearity

$\frac{{\text d}y}{{\text d}x} = \cos \left( 2y-x \right) +1 , \qquad y(0) = \pi /4 .$
This problem can be solve at least implicitly by substitution: z = 2y - x. Then for new dependent variable, we get a separable equation
$z' = 2y' -1 = 2\,\cos z -1 , z(0) = \pi /2 .$
Separating variables and integrating, we obtain
$\frac{{\text d}z}{2\,\cos z - 1} = - {\text d}x \qquad \Longrightarrow \qquad \frac{2}{\sqrt{3}}\,\mbox{arctanh} \left[ \sqrt{3}\,\tan \frac{z}{2} \right] = x .$
Integrate[1/(2*Cos[z] - 1), z]
(2 ArcTanh[Sqrt Tan[z/2]])/Sqrt
Now we can expand the solution into power series:
$y(x) = \frac{\pi}{2} - x + \frac{4}{3!}\,x^3 + \cdots . \qquad x\ne -1.$
Now we use the modified decomposition method, assuming that the solution is represented by a convergent series
$y(x) = \sum_{n\ge 0} c_n x^n = \frac{\pi}{2} + \sum_{n\ge 1} c_n x^n .$
With the initial term c0 = π/2, we calculated the Adomian polynomials:
\begin{align*} A_0 &= \sin \left( x - u_0 \right) = \sin \left( x - \frac{\pi}{2} \right) = - \cos x = -1 + \frac{x^2}{2} - \frac{x^4}{4!} + \cdots , \\ A_1 &= \end{align*}
■

Example: Consider the initial value problem

$y' = \frac{3}{8}\, y - \frac{9t}{8\,y} , \qquad y(0) = 2,$
that has the explicit solution
$y(t) = \sqrt{4+3t} = 2 + \frac{3}{4}\, t + \left( \frac{3\,t}{8} \right)^2 + \left( \frac{3\,t}{8} \right)^3 + \frac{5}{4} \left( \frac{3\,t}{8} \right)^4 + \frac{7}{4} \left( \frac{3\,t}{8} \right)^5 + \frac{21}{8} \left( \frac{3\,t}{8} \right)^6 + \frac{33}{8} \left( \frac{3\,t}{8} \right)^7 + \frac{429}{8^2} \left( \frac{3\,t}{8} \right)^8 + \cdots .$
Series[Sqrt[4 + 3*t], {t, 0, 8}]
$$y(t) = \sqrt{4+3t} .$$    ■

1. Adomian, G., Rach, R., Transformation of series, Applied Mathematics Letters Volume 4, Issue 4, 1991, Pages 69-71; https://doi.org/10.1016/0893-9659(91)90058-4
2. Adomian, G., Rach, R., Nonlinear transformation of series—part II, Computers & Mathematics with Applications Volume 23, Issue 10, May 1992, Pages 79-83; https://doi.org/10.1016/0898-1221(92)90058-P
3. Andrianov, I.V., Olevskii, V.I., Tokarzewski, S., A modified Adomian's decomposition method, Journal of Applied Mathematics and Mechanics, 1998, Volume 62, Issue 2, Pages 309--314; https://doi.org/10.1016/S0021-8928(98)00040-9
4. Duan, J.-S. and Rach, R., The degenerate form of the Adomian polynomials in the power series method for nonlinear ordinary differential equations, Journal of Mathematics and System Science, Volume 5, Pages 411--428, doi: 10.17265/2159-5291/2015.10.003
5. Duan, J.-S., Rach, R., Wazwaz, A.-M., Solution of the model of beam-type micro- and nano-scale electrostatic actuators by a new modified Adomian decomposition method for nonlinear boundary value problems, International Journal of Non-Linear Mechanics, 2013, Volume 49, 2013, Pages 159--169.
6. Rach, R., Adomian, G., Meyers, R., A modified decomposition method, Computers & Mathematics with Applications, 1992, Volume 23, Issue 1, January 1992, Pages 17-23; https://doi.org/10.1016/0898-1221(92)90076-T