**3. Basic idea of Homotopy perturbation method**

First, we outline the general procedure of the homotopy perturbation method developed and advanced by He. We consider the differential Eq. [2]

$$A(u) - f\left(r\right) = 0, r \in \Omega\tag{4}$$

$$B\left(u,\frac{\partial u}{\partial \mathbf{x}}\right) = \mathbf{0}, r \in \Gamma \tag{5}$$

where *A* is a general differential operator, linear or nonlinear, *f r*ð Þ is a known analytic function, *B* is a boundary operator and Γ is the boundary of the domain Ω. The operator *A* can be generally divided into two operators, *L* and *N*, where *L* is linear and *N* is a nonlinear operator. Eq. (4) can be written as

$$L(u) + N(u) - f\left(r\right) = 0\tag{6}$$

Using the homotopy technique, we can construct a homotopy [1,2]

$$v(r, p): \Omega \times [0, 1] \to \mathbb{R} \text{ which satisfies the relation}$$

$$H(v, p) = (\mathbb{1} - p)[L(v) - L(u\_0)] + p[A(v) - f(r)] = \mathbf{0}, r \in \Omega \tag{7}$$

Here *p* ∈½ � 0, 1 is called the homotopy parameter and *u*<sup>0</sup> is an initial approximation for the solution of Eq. (4), which satisfies the boundary conditions. Clearly, from Eq. (7), we have

$$H(\upsilon, \mathbf{0}) = L(\upsilon) - L(\mu\_0) \tag{8}$$

$$H(v, \mathbf{1}) = A(v) - f\left(r\right) \tag{9}$$

We assume that the solution of Eq. (7) can be expressed as a series in *p* as follows:

$$\nu = \nu\_0 + p\nu\_0 + p\_2\nu\_2 + p\_3\nu\_3 + \dotsb \tag{10}$$

On setting *p* ¼ 1, we obtain the approximate solution of Eq. (10) as

$$\mu = \lim\_{p \to 1} \nu = \nu\_0 + \nu\_0 + \nu\_2 + \nu\_3 + \dotsb \tag{11}$$
