**3. Solution Method and Details**

In order to solve Equations (14–18) under the boundary conditions (19, 20), we use the Homotopy Analysis Method (HAM) with the following procedure. The solutions with the auxiliary parameters ¯*h* adjust and control the convergence of the solutions.

The initial guesses are selected as follows:

$$f\_0(\eta) = \begin{pmatrix} 1 - \varepsilon^{-\eta} \end{pmatrix}, \; \theta\_{1,0}(\eta) = \varepsilon^{-\eta}, \; \theta\_{2,0}(\eta) = \eta \varepsilon^{-\eta}, \; \chi\_{1,0}(\eta) = \varepsilon^{-\eta}, \; \chi\_{2,0}(\eta) \tag{26}$$

The linear operators are taken as *Lf* , *Lθ*<sup>1</sup> , *Lθ*<sup>2</sup> , *Lχ*<sup>1</sup> , *Lχ*<sup>2</sup>

$$\begin{array}{c} L\_f(f) = f''' - f', \ L\_{\theta\_1}(\theta\_1) = \theta\_1'' - \theta\_1, \ L\_{\theta\_2}(\theta\_2) = \theta\_2'' - \theta\_2\\ L\_{\chi\_1}(\chi\_1) = \chi\_1'' - \chi\_1, \ L\_{\chi\_2}(\chi\_2) = \chi\_2'' - \chi\_2 \end{array} \tag{27}$$

which have the following properties:

$$\begin{array}{c} L\_f(\mathfrak{c}\_1 + \mathfrak{c}\_2 \mathfrak{e}^{-\eta} + \mathfrak{c}\_3 \mathfrak{e}^{\eta}) = 0, \ L\_{\mathfrak{\theta}\_1}(\mathfrak{c}\_4 \mathfrak{e}^{\eta} + \mathfrak{c}\_5 \mathfrak{e}^{-\eta}) = 0\\ L\_{\mathfrak{\theta}\_2}(\mathfrak{c}\_6 \mathfrak{e}^{\eta} + \mathfrak{c}\_7 \mathfrak{e}^{-\eta}) = 0, \ L\_{\chi\_1}(\mathfrak{c}\_8 \mathfrak{e}^{-\eta} + \mathfrak{c}\_9 \mathfrak{e}^{\eta}) = 0, \ L\_{\chi\_1}(\mathfrak{c}\_{10} \mathfrak{e}^{-\eta} + \mathfrak{c}\_{11} \mathfrak{e}^{\eta}) = 0 \end{array} \tag{28}$$

where *ci*(*i* = 1 − 11) are the constants in general solution: The resultant non-linear operatives *Nf* , *Nθ*<sup>1</sup> , *Nθ*<sup>2</sup> , *Nχ*<sup>1</sup> , *Nχ*<sup>2</sup> are given as

$$N\_f[f(\eta; p), \theta\_1(\eta; p)] = \left(\varepsilon\_2 \mathcal{W} \varepsilon + \frac{k\_{\text{Brownian}}}{k\_f Pr\_f \varepsilon\_2}\right) \frac{\partial^3 f(\eta; p)}{\partial \eta^3} - \varepsilon \varepsilon \frac{2\mathcal{\beta}}{(\eta + a)^4} \theta\_1(\eta; p) - \left(\frac{\partial f(\eta; p)}{\partial \eta}\right)^2 + f(\eta; p) \frac{\partial^2 f(\eta; p)}{\partial \eta^2} \tag{29}$$

*Nθ*<sup>1</sup> [ *f*(*η*; *p*), *θ*1(*η*; *p*), *θ*2(*η*; *p*)] = *ε*<sup>2</sup> *kn f k f* 1 *Pr* -*∂*2*θ*1(*η*;*p*) *∂η*<sup>2</sup> + <sup>2</sup>*θ*2(*η*; *<sup>p</sup>*) + *f*(*η*; *p*) *∂θ*1(*η*;*p*) *∂η* <sup>+</sup> *<sup>ε</sup>*<sup>2</sup> <sup>1</sup> *Pr* <sup>2</sup>*<sup>λ</sup>* (*η*+*α*) <sup>4</sup> (*f*(*η*; *p*)*θ*1(*η*; *p*) − *εf*(*η*; *p*))− *δe* - (*f*(*η*; *<sup>p</sup>*))<sup>2</sup> *<sup>∂</sup>*2*θ*1(*η*;*p*) *∂η*<sup>2</sup> <sup>+</sup> *<sup>f</sup>*(*η*; *<sup>p</sup>*) *<sup>∂</sup> <sup>f</sup>*(*η*;*p*) *∂η ∂θ*1(*η*;*p*) *∂η* (30) *Nθ*<sup>2</sup> [ *f*(*η*; *p*), *θ*1(*η*; *p*), *θ*2(*η*; *p*)] = *ε*<sup>2</sup> *kn f k f* 1 *Pr ∂*2*θ*2(*η*;*p*) *∂η*<sup>2</sup> + *λ Pr* - *ε*2*We* + *ε*<sup>1</sup> *kBrownian k <sup>f</sup> Prf ε*<sup>2</sup> *f*(*η*; *p*) -*∂ f*(*η*;*p*) *∂η* <sup>2</sup> *θ*1(*η*; *p*) + *ε*<sup>2</sup> 2*λ Pr* <sup>1</sup> (*η*+*α*) <sup>3</sup> *f*(*η*; *p*)*θ*2(*η*; *p*) −*ε*<sup>2</sup> *λβ Pr* 4 (*η*+*α*) <sup>5</sup> *f*(*η*; *p*) + <sup>2</sup> (*η*+*α*) 4 *∂ f*(*η*;*p*) *∂η δe* ⎛ ⎝ (*f*(*η*; *<sup>p</sup>*))<sup>2</sup> *<sup>∂</sup>*2*θ*2(*η*;*p*) *∂η*<sup>2</sup> <sup>−</sup> <sup>3</sup> *<sup>f</sup>*(*η*; *<sup>p</sup>*) *<sup>∂</sup> <sup>f</sup>*(*η*;*p*) *∂η ∂θ*2(*η*;*p*) *∂η* <sup>−</sup>*f*(*η*; *<sup>p</sup>*) *<sup>∂</sup>*<sup>2</sup> *<sup>f</sup>*(*η*;*p*) *∂η*<sup>2</sup> *<sup>θ</sup>*2(*η*; *<sup>p</sup>*) + <sup>4</sup> -*∂ f*(*η*;*p*) *∂η* <sup>2</sup> *θ*2(*η*; *p*) ⎞ ⎠ (31)

$$\begin{split} N\_{\chi\_1}[f(\eta;p), \chi\_1(\eta;p), \chi\_2(\eta;p)] &= (1-\phi)^{2.5} \frac{1}{\text{Sc}} \left( \frac{\partial^2 \chi\_1(\eta;p)}{\partial \eta^2} + 2\chi\_2(\eta;p) \right) + \\ f(\eta;p) \frac{\partial \chi\_1(\eta;p)}{\partial \eta} - \sigma \chi\_1(\eta;p) \end{split} \tag{32}$$

$$\begin{split} N\_{\chi\_2}[f(\eta;p),\chi\_2(\eta;p)] &= (1-\phi)^{2.5} \frac{1}{\mathbb{S}\varepsilon} \frac{\partial^2 \chi\_1(\eta;p)}{\partial \eta^2} + f(\eta;p) \frac{\partial \chi\_2(\eta;p)}{\partial \eta} \\ &- 2 \frac{\partial f(\eta;p)}{\partial \eta} \chi\_2(\eta;p) - \sigma \chi\_2(\eta;p) \end{split} \tag{33}$$

The basic idea of the HAM is described in [1–7]; the zeroth-order problems from Equations (14)–(18) are

$$(1-p)L\_f[f(\eta;p)-f\_0(\eta)] = p\hbar\_f N\_f[f(\eta;p), \theta\_1(\eta;p)] \tag{34}$$

$$\left[ (1-p)L\_{\theta\_1} \left[ \theta\_1(\eta; p) - \theta\_{1,0}(\eta) \right] \right] = p\hbar\_{\theta\_1} \mathcal{N}\_{\theta\_1} \left[ f(\eta; p), \theta\_1(\eta; p), \theta\_2(\eta; p) \right] \tag{35}$$

$$\mathbb{P}\left(1-p\right)L\_{\theta\_2}\left[\theta\_2(\eta;p)-\theta\_{2,\,0}(\eta)\right]=p\hbar\_{\theta\_2}\mathcal{N}\_{\theta\_2}\left[f(\eta;p),\,\theta\_1(\eta;p),\,\theta\_2(\eta;p)\right] \tag{36}$$

$$[(1-p)L\_{\chi\_1}[\chi\_1(\eta;p)-\chi\_{1,0}(\eta)] = p\hbar\_{\chi\_1}N\_{\chi\_1}[f(\eta;p),\chi\_1(\eta;p),\ \chi\_2(\eta;p)] \tag{37}$$

$$(1 - p)L\_{\chi1}[\chi\_2(\eta; p) - \chi\_{2, \, 0}(\eta)] = p\hbar\_{\chi2}N\_{\chi2}[f(\eta; p), \chi\_2(\eta; p)]\tag{38}$$

The equivalent boundary conditions are

$$\begin{array}{ll} f(\eta;p)|\_{\eta=0} = 0, & \frac{\partial f(\eta;p)}{\partial \eta} \Big|\_{\eta=0} = 1, & \frac{\partial f(\eta;p)}{\partial \eta} \Big|\_{\eta \to \infty} = 0 \\ & & \theta\_1(\eta;p)|\_{\eta=0} = 1, & \theta\_1(\eta;p)|\_{\eta \to \infty} = 0 \\ & & \theta\_2(\eta;p)|\_{\eta=0} = 0, & \theta\_2(\eta;p)|\_{\eta \to \infty} = 0 \\ & & \chi\_1(\eta;p)|\_{\eta=0} = 1, & \chi\_1(\eta;p)|\_{\eta \to \infty} = 0 \\ & & \chi\_2(\eta;p)|\_{\eta=0} = 0, & \chi\_2(\eta;p)|\_{\eta \to \infty} = 0 \end{array} \tag{39}$$

where *p* ∈ [0, 1] is the imbedding parameter; *h*¯ *<sup>f</sup>* , *h*¯ *<sup>θ</sup>*<sup>1</sup> , *h*¯ *<sup>θ</sup>*<sup>2</sup> , *h*¯ *<sup>χ</sup>*<sup>1</sup> , *h*¯ *<sup>χ</sup>*<sup>2</sup> are used to control the convergence of the solution. When *p* = 0 and *p* = 1, we have

$$f(\eta;1) = f(\eta), \; \theta\_1(\eta;1) = \theta\_1(\eta), \; \theta\_2(\eta;1) = \theta\_2(\eta), \; \chi\_1(\eta;1) = \chi\_1(\eta), \; \chi\_2(\eta;1) = \chi\_2(\eta) \tag{40}$$

Expanding *f*(*η*; *p*), *θ*1(*η*; *p*), *θ*2(*η*; *p*), *χ*1(*η*; *p*), *χ*2(*η*; *p*) in Taylor's series about *p* = 0

$$\begin{aligned} f(\eta; p) &= f\_0(\eta) + \sum\_{\substack{\mathfrak{m} = 1 \\ \mathfrak{m} = 1 \\ \theta\_1(\eta; p) = \theta\_{1,0}(\eta) + \sum\_{\substack{\mathfrak{m} = 1 \\ \mathfrak{m} = 1}}^{\infty} \theta\_{1,\mathfrak{m}}(\eta) p^{\mathfrak{m}} \\ \theta\_2(\eta; p) = \theta\_{2,0}(\eta) + \sum\_{\substack{\mathfrak{m} = 1 \\ \mathfrak{m} = 1 \\ \mathfrak{m} = 1}}^{\infty} \theta\_{2,\mathfrak{m}}(\eta) p^{\mathfrak{m}} \\ \chi\_1(\eta; p) &= \chi\_{1,0}(\eta) + \sum\_{\substack{\mathfrak{m} = 1 \\ \mathfrak{m} = 1 \\ \mathfrak{m} = 1}}^{\infty} \chi\_{1,\mathfrak{m}}(\eta) p^{\mathfrak{m}} \\ \chi\_2(\eta; p) &= \chi\_{2,0}(\eta) + \sum\_{\mathfrak{m} = 1}^{\infty} \chi\_{2,\mathfrak{m}}(\eta) p^{\mathfrak{m}} \end{aligned} \tag{41}$$

where

$$f\_{\mathfrak{m}}(\eta) = \frac{1}{\stackrel{\mathfrak{m}}{m!}} \frac{\partial f(\eta; p)}{\partial \eta} \Big|\_{p=0'} \theta\_{1, \mathfrak{m}}(\eta) = \frac{1}{\stackrel{\mathfrak{m}!}{m!}} \frac{\partial \theta\_1(\eta; p)}{\partial \eta} \Big|\_{p=0}$$

$$\theta\_{2, \mathfrak{m}}(\eta) = \frac{1}{\stackrel{\mathfrak{m}!}{m!}} \frac{\partial \theta\_2(\eta; p)}{\partial \eta} \Big|\_{p=0'} \chi\_{1, \mathfrak{m}}(\eta) = \frac{1}{\stackrel{\mathfrak{m}!}{m!}} \frac{\partial \chi\_1(\eta; p)}{\partial \eta} \Big|\_{p=0'} \chi\_{2, \mathfrak{m}}(\eta) = \frac{1}{\stackrel{\mathfrak{m}!}{m!}} \frac{\partial \chi\_2(\eta; p)}{\partial \eta} \Big|\_{p=0} \tag{42}$$

The secondary constraints *h*¯ *<sup>f</sup>* , *h*¯ *<sup>θ</sup>*<sup>1</sup> , *h*¯ *<sup>θ</sup>*<sup>2</sup> , *h*¯ *<sup>χ</sup>*<sup>1</sup> , *h*¯ *<sup>χ</sup>*<sup>2</sup> are chosen in such a way that the series (40) converges at *p* = 1; switching *p* = 1 in (40), we obtain

$$\begin{aligned} f(\eta) &= f\_0(\eta) + \sum\_{m=1}^{\infty} f\_m(\eta) \\ \theta\_1(\eta) &= \theta\_{1,0}(\eta) + \sum\_{\substack{m=1 \\ m=1}}^{\infty} \theta\_{1,m}(\eta) \\ \theta\_2(\eta) &= \theta\_{2,0}(\eta) + \sum\_{m=1}^{\infty} \theta\_{2,m}(\eta) \\ \chi\_1(\eta) &= \chi\_{1,0}(\eta) + \sum\_{\substack{m=1 \\ m=1 \\ m=1}}^{\infty} \chi\_{1,m}(\eta) \\ \chi\_2(\eta) &= \chi\_{2,0}(\eta) + \sum\_{m=1}^{\infty} \chi\_{2,m}(\eta) \end{aligned} \tag{43}$$

The *mth*-*order* problem satisfies the following:

$$\begin{aligned} L\_f[f\_m(\eta) - \chi\_m f\_{m-1}(\eta)] &= \hbar\_f R\_m^f(\eta) \\ L\_{\theta\_1}[\theta\_{1,m}(\eta) - \chi\_m \theta\_{1,m-1}(\eta)] &= \hbar\_{\theta\_1} \mathbb{R}\_m^{\theta\_1}(\eta) \\ L\_{\theta\_2}[\theta\_{2,m}(\eta) - \chi\_m \theta\_{2,m-1}(\eta)] &= \hbar\_{\theta\_2} R\_m^{\theta\_2}(\eta) \\ L\_{\chi\_1}[\chi\_{1,m}(\eta) - \chi\_m \chi\_{1,m-1}(\eta)] &= \hbar\_{\chi\_1} R\_m^{\chi\_1}(\eta) \\ L\_{\chi\_2}[\chi\_{2,m}(\eta) - \chi\_m \chi\_{2,m-1}(\eta)] &= \hbar\_{\chi\_2} R\_m^{\chi\_2}(\eta) \end{aligned} \tag{44}$$

The corresponding boundary conditions are as follows:

$$\begin{aligned} f\_m(0) &= f\_m'(0) = \theta\_{1,m}'(0) = \theta\_{2,m}'(0) = \chi\_{1,m}(0) = \chi\_{2,m}(0) = 0 \\ f\_m'(\infty) &= \theta\_{1,m}(\infty) = \theta\_{2,m}(\infty) = \chi\_{1,m}(\infty) = \chi\_{2,m}(\infty) = 0 \end{aligned} \tag{45}$$

Here

$$\begin{split} R\_{m}^{f}(\eta) &= \left(\varepsilon\_{2}\mathcal{W}\varepsilon + \frac{k\_{\text{Brownian}}}{k\_{f}Pr\_{f}\varepsilon\_{2}}\right)f\_{m-1}^{\prime\prime} - \varepsilon\_{2}\frac{2\delta}{(\eta+a)^{4}}\theta\_{1,m-1} - \frac{\partial^{3}f(\eta;p)}{\partial\eta^{3}} -\\ &\varepsilon\_{2}\frac{2\delta}{(\eta+a)^{4}}\theta\_{1}(\eta;p) - \sum\_{k=0}^{m-1}f\_{m-1}^{\prime}f\_{k}^{\prime} + \sum\_{k=0}^{m-1}f\_{m-1-k}f\_{k}^{\prime\prime} \end{split} \tag{46}$$

*Rθ*1 *<sup>m</sup>* (*η*) = *ε*<sup>2</sup> *kn f k f* 1 *Pr* - *θ* 1, *<sup>m</sup>*−<sup>1</sup> <sup>+</sup> <sup>2</sup>*θ*2, *<sup>m</sup>*−<sup>1</sup> + *ε*<sup>2</sup> <sup>1</sup> *Pr* <sup>2</sup>*<sup>λ</sup>* (*η*+*α*) 4 *m*−1 ∑ *k*=0 *fm*−1−*kθ*1, *<sup>k</sup>* − *εfm*−<sup>1</sup> <sup>+</sup> *<sup>m</sup>*−<sup>1</sup> ∑ *k*=0 *fm*−1−*kθ* 1, *<sup>k</sup>* − *δ<sup>e</sup> m*−1 ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *fk*−*iθ* 1,*<sup>i</sup>* <sup>+</sup> *<sup>m</sup>*−<sup>1</sup> ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *f k*−*i θ*1, *<sup>i</sup>* (47) *Rθ*<sup>2</sup> *<sup>m</sup>* (*η*) = *ε*<sup>2</sup> *kn f k f* 1 *Pr θ* 2,*m*−<sup>1</sup> <sup>+</sup> *<sup>λ</sup> Pr* - *ε*2*We* + *ε*<sup>1</sup> *kBrownian k <sup>f</sup> Prf ε*<sup>2</sup> *<sup>m</sup>*−<sup>1</sup> ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *f k i* ∑ *p*=0 *f <sup>i</sup>*−*pθ*1,*<sup>p</sup>* −*ε*<sup>2</sup> *λβ Pr* 4 (*η*+*α*) <sup>5</sup> *fm*−<sup>1</sup> <sup>+</sup> <sup>2</sup> (*η*+*α*) <sup>4</sup> *f m*−1 + *ε*<sup>2</sup> 2*λβ Pr* <sup>1</sup> (*η*+*α*) 3 *m*−1 ∑ *k*=0 *fm*−1−*kθ*2,*<sup>k</sup> δe* ⎡ ⎢ ⎢ ⎣ *m*−1 ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *fk*−*iθ* 2,*<sup>i</sup>* − 3 *m*−1 ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *f k*−*i θ*2,*i*− *m*−1 ∑ *k*=0 *fm*−1−*<sup>k</sup> k* ∑ *i*=0 *f k*−*i θ*2,*<sup>i</sup>* + 4 *m*−1 ∑ *k*=0 *f m*−1−*k k* ∑ *i*=0 *f k*−*i θ*2,*<sup>i</sup>* ⎤ ⎥ ⎥ ⎦ (48) *Rχ*<sup>1</sup> 2.5 1 - *χ* + *m*−1 ∑

$$R\_{\mathfrak{m}}^{\lambda\_1}(\eta) = (1 - \phi)^{2.5} \frac{1}{Sc} \left( \chi\_{1, \mathfrak{m}-1}^{\eta} + 2\chi\_{2, \mathfrak{m}-1} \right) + \sum\_{k=0}^{\cdot} f\_{\mathfrak{m}-1-k} \chi\_{1,k}^{\prime} - \sigma \chi\_{1, \mathfrak{m}-1} \tag{49}$$
 
$$\chi\_{1, \mathfrak{m}-1, \mathfrak{m}}^{\prime} = \dots \quad \chi\_{2.5}^{\prime} \stackrel{\mathfrak{m}-1}{\longrightarrow} \dots \quad \stackrel{\mathfrak{m}-1}{\longrightarrow} \dots \quad \stackrel{\mathfrak{m}-1}{\longrightarrow} \chi\_{1, \mathfrak{m}-1, \mathfrak{m}-1} \tag{50}$$

$$R\_m^{X\_2}(\eta) = (1 - \phi)^{2.5} \frac{1}{Sc} \chi\_{2, m - 1}^{''} + \sum\_{k = 0}^{m - 1} f\_{m - 1 - k} \chi\_{2, k}^{'} - \sum\_{k = 0}^{m - 1} f\_{m - 1 - k}^{'} \chi\_{2, k} - \sigma \chi\_{2, m - 1} \tag{50}$$

where

$$\chi\_m = \begin{cases} \begin{array}{c} 0, \text{ } if \ p \le 1 \\ 1, \text{ } if \ p > 1 \end{array} \tag{51}$$

### *Validation and Comparison*

Table 1 shows the physical properties of nanofluid. A comparison of the validation of the results for the velocity, temperature, and concentration fields using the Homotopy Analysis Method and a numerical (ND-solve) method are shown in Tables 2–4 and in Figures 2–4. From these tables and figures, it can be observed that the results of both methods are in good agreement.



**Table 2.** Comparison table for HAM solution and numerical method for velocity field and their results.


**Table 3.** Comparison table for HAM solution and numerical method for temperature field and their results.

**Table 4.** Comparison table for HAM solution and numerical method for concentration field and their results.


**Figure 2.** Comparison graph for velocity profile.

**Figure 3.** Comparison graph for temperature profile.

**Figure 4.** Comparison graph for concentration profile.
