**2. Problem Description**

A non-linear problem, concerning the transport and heat-transfer characteristics of non-Newtonian (Ag-TiO2/H2O) hybrid nanofluid in an endoscope due to ciliary metachronical rhythm is investigated here. A cylindrical coordinate system is used with (*R*, θ,*Z*) as position coordinates of fluid particles. Non-Newtonian behavior of the flow is considered with Ostwald-de-Waele power law model expressed as [22,25]:

$$\tau = -k \left\{ \sqrt{\left| \frac{1}{2} \Delta : \Delta \right|}^{n-1} \right\} \Lambda. \tag{1}$$

where,

$$\frac{1}{2}\Delta : \Delta = 2\left(\left(\frac{\partial \mathcal{U}}{\partial \mathcal{R}}\right)^2 + \left(\frac{\mathcal{U}}{\mathcal{R}}\right)^2 + \left(\frac{\partial \mathcal{W}}{\partial \mathcal{Z}}\right)^2\right) + \left(\frac{\partial \mathcal{W}}{\partial \mathcal{R}} + \frac{\partial \mathcal{U}}{\partial \mathcal{Z}}\right)^2. \tag{2}$$

In above expression *k* and n respectively represents flow consistency index and power law index. In this study, non-Newtonian shear thickening fluid is considered with *n* = 2. Moreover, a constant magnetic field having intensity H0 is taken in radial direction causes an induced magnetic field H- (Hr(r, z), 0, Hz(r, z)) and, thereby, total magnetic field vector is H(H0 + Hr(r, z), 0, Hz(r, z)).

The governing system of equations in an unsteady form is [25,29,30]:

$$\frac{\partial \mathbf{U}}{\partial \mathbf{R}} + \frac{\mathbf{U}}{\mathbf{R}} + \frac{\partial \mathbf{W}}{\partial \mathbf{Z}} = \mathbf{0} \tag{3}$$

$$\frac{\partial \mathcal{U}}{\partial t} + \mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathbf{K}} + \mathcal{W}\frac{\partial \mathcal{U}}{\partial \mathcal{Z}} = \frac{-1}{\rho\_{\text{fluid}}}\frac{\partial \mathcal{P}}{\partial \mathcal{R}} - \frac{1}{\rho\_{\text{fluid}}} \left( \frac{\partial \left( 2\text{k}\text{R}\mathbf{o}\frac{\partial \mathcal{U}}{\partial \mathbf{K}} \right)}{\partial \mathbf{K}} + \frac{\partial \left( \text{k}\mathbb{1}\left( \frac{\partial \mathcal{U}}{\partial \mathbf{Z}} + \frac{\partial \mathcal{U}}{\partial \mathbf{K}} \right) \right)}{\partial \mathbf{K}} \right) - \frac{\partial}{2\rho\_{\text{fluid}}} \left( \frac{\partial \mathcal{U}}{\partial \mathbf{K}} \right) + \mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathcal{R}} + \mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathcal{R}} + \mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathcal{R}}$$
 
$$\frac{\partial}{\partial \mathbf{m}} \left( \frac{\partial \mathcal{U}\_{\mathbf{r}}}{\partial \mathbf{t}} + (\mathbf{H}\_{0} + \mathbf{H}\_{\mathbf{r}}) \frac{\partial \mathcal{H}\_{\mathbf{r}}}{\partial \mathbf{K}} + \mathbf{H}\_{\mathbf{r}}\frac{\partial \mathcal{U}\_{\mathbf{r}}}{\partial \mathcal{Z}} \right),$$

$$
\frac{\partial \mathbf{W}}{\partial t} + \mathbf{U} \frac{\partial \mathbf{W}}{\partial \mathbf{R}} + \mathbf{W} \frac{\partial \mathbf{W}}{\partial \mathbf{Z}} = \frac{-1}{\rho\_{\text{fml}}} \frac{\partial \mathbf{P}}{\partial \mathbf{Z}} - \frac{1}{\rho\_{\text{fml}}} \left( \frac{\partial \left( \mathbf{k} \otimes \mathbf{R} \left( \frac{\partial \mathbf{J}}{\partial \mathbf{Z}} + \frac{\partial \mathbf{W}}{\partial \mathbf{R}} \right) \right)}{\partial \mathbf{R}} + \frac{\partial \left( \mathbf{k} \otimes \mathbf{r} \frac{\partial \mathbf{V}}{\partial \mathbf{Z}} \right)}{\partial \mathbf{R}} \right) - \frac{\partial}{2\rho\_{\text{fml}}} \left( \frac{\partial \mathbf{H}}{\partial \mathbf{Z}} \right) + \mathbf{M} \tag{5}
$$

$$
\frac{\frac{\partial}{\partial \mathbf{m}} \left( \frac{\partial \mathbf{H}\_{\mathbf{Z}}}{\partial \mathbf{R}} + (\mathbf{H}\_{0} + \mathbf{H}\_{\mathbf{r}}) \frac{\partial \mathbf{H}\_{\mathbf{Z}}}{\partial \mathbf{R}} + \mathbf{H}\_{\mathbf{Z}} \frac{\partial \mathbf{H}\_{\mathbf{Z}}}{\partial \mathbf{Z}} \right) + \frac{(\rho \partial)\_{\text{fml}}}{\rho\_{\text{fml}}} \mathbf{g} (\mathbf{T} - \mathbf{T}\_{0}),
$$

$$\frac{\partial \mathbf{T}}{\partial \mathbf{t}} + \mathbf{U} \frac{\partial \mathbf{T}}{\partial \mathbf{R}} + \mathbf{W} \frac{\partial \mathbf{T}}{\partial \mathbf{Z}} = \frac{\kappa\_{\text{hnf}}}{\left(\rho \mathbf{c}\_{\text{P}}\right)\_{\text{hnf}}} \left(\frac{\partial^{2} \mathbf{T}}{\partial \mathbf{R}^{2}} + \frac{1}{\mathbf{R}} \frac{\partial \mathbf{T}}{\partial \mathbf{R}} + \frac{\partial^{2} \mathbf{T}}{\partial \mathbf{Z}^{2}}\right) + \frac{\mathbf{Q}\_{0}}{\left(\rho \mathbf{c}\_{\text{P}}\right)\_{\text{hnf}}},\tag{6}$$

$$\begin{split} \frac{1}{\sqrt{\mu}} \frac{\partial \overline{\mathbf{L}}}{\partial \overline{\mathbf{L}}} &= \mathbf{U} \left( -\frac{\partial \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{R}}} + \frac{\partial \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{Z}}} \right) + \mathbf{H}\_{\mathbf{z}} \frac{\partial \mathbf{U}}{\partial \overline{\mathbf{Z}}} - \left( \mathbf{H}\_{0} + \mathbf{H}\_{\mathbf{r}} \right) \frac{\partial \mathbf{V}}{\partial \overline{\mathbf{Z}}} - 2\mathcal{W} \frac{\partial \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{Z}}} + \\\\ \frac{1}{\sigma \overline{\theta}} \left( \frac{\partial^{2} \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{R}}^{2}} + \frac{1}{\overline{\mathbf{R}}} \frac{\partial \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{R}}} + \frac{\partial^{2} \mathbf{H}\_{\mathbf{r}}}{\partial \overline{\mathbf{Z}}^{2}} \right). \end{split} \tag{7}$$

$$\begin{split} \frac{-1}{\rho} \frac{\partial \mathbb{E}}{\partial \mathbb{R}} &= \quad (\mathbb{H}\_0 + \mathbb{H}\_1) \frac{\partial \mathbb{W}}{\partial \mathbb{R}} + \mathbb{W} \Big( \frac{\partial \mathbb{H}\_1}{\partial \mathbb{R}} - \frac{\partial \mathbb{H}\_1}{\partial \mathbb{Z}} \Big) - 2 \mathbb{U} \frac{\partial \mathbb{H}\_1}{\partial \mathbb{R}} - \mathbb{H}\_2 \frac{\partial \mathbb{U}}{\partial \mathbb{R}} + \\\\ &\frac{1}{\sigma \overline{\rho}} \Big( \frac{\partial^2 \mathbb{H}\_x}{\partial \mathbb{R}^2} + \frac{1}{\mathbb{R}} \frac{\partial \mathbb{H}\_x}{\partial \mathbb{R}} + \frac{\partial^2 \mathbb{H}\_x}{\partial \mathbb{Z}^2} \Big). \end{split} \tag{8}$$

where,

$$\mathcal{O} = \sqrt{\left| 2 \overline{\left[ \left( \frac{\partial \mathcal{U}}{\partial \mathcal{R}} \right)^2 + \left( \frac{\mathcal{U}}{\mathcal{R}} \right)^2 + \left( \frac{\partial \mathcal{W}}{\partial \mathcal{Z}} \right)^2 \right] + \left( \frac{\partial \mathcal{U}}{\partial \mathcal{Z}} + \frac{\partial \mathcal{W}}{\partial \mathcal{R}} \right)^2 \right|^{n-1}} $$

Wave shapes in the laboratory frame for envelope of cilia tips according to Figure 1 can be expressed as:

**Figure 1.** Geometry of the physical problem.

$$R\_1 = a\_{1\prime} \tag{9}$$

$$R\_2 = f(Z, t) = \left[ a\_2 + b \cos\left(\frac{2\pi}{\lambda}(Z - ct)\right) \right] \tag{10}$$

where, α<sup>1</sup> and α<sup>2</sup> represents radii of internal and external cylindrical tubes, accordingly. Considering the motion of cilia in an elliptical path, the vertical position of cilia tips is expressed as:

$$Z = \mathbf{g}(Z, Z\_0, t) = \left[ Z\_0 + \alpha \text{bSim}(\frac{2\pi}{\lambda}(Z - \text{ct})) \right],\tag{11}$$

Since the velocities of the fluid layers are similar to those of the cilia tips under the no-slip wall conditions, the vertical and horizontal velocities are:

$$\begin{split} \mathbf{W} &= \left. \frac{\partial \mathbf{Z}}{\partial t} \right|\_{\mathbf{Z}\_0} = \frac{\partial \boldsymbol{\xi}}{\partial t} + \frac{\partial \boldsymbol{\xi}}{\partial \mathbf{Z}} \frac{\partial \mathbf{Z}}{\partial t}, \\ \mathbf{U} &= \left. \frac{\partial \mathbf{R}}{\partial t} \right|\_{\mathbf{Z}\_0} = \frac{\partial f}{\partial t} + \frac{\partial f}{\partial \mathbf{Z}} \frac{\partial \mathbf{Z}}{\partial t}. \end{split} \tag{12}$$

With the help of Equations (10) and (11), Equation (12) becomes:

$$\text{At R}=\text{R}\_2\text{ W}=\frac{\frac{-2\pi}{\lambda}\text{bar}\text{Cov}\left(\frac{2\pi}{\lambda}(Z-\text{ct})\right)}{1-\frac{2\pi}{\lambda}\text{bar}\text{Cov}\left(\frac{2\pi}{\lambda}(Z-\text{ct})\right)},\text{ U}=\frac{\frac{2\pi}{\lambda}\text{bar}\text{Cov}\left(\frac{2\pi}{\lambda}(Z-\text{ct})\right)}{1-\frac{2\pi}{\lambda}\text{bar}\text{Cov}\left(\frac{2\pi}{\lambda}(Z-\text{ct})\right)}.\tag{13}$$

The associated boundary conditions are defined as:

$$\mathbf{W} = \mathbf{0}, \text{ at } \mathbf{R} = \mathbf{R}\_1, \mathbf{W} = \frac{\frac{-2\pi}{\lambda} \mathbf{b} \alpha \mathbf{c} \mathrm{Cos} \left(\frac{2\pi}{\lambda} (\mathbf{Z} - \mathbf{c}\mathbf{t})\right)}{1 - \frac{2\pi}{\lambda} \mathbf{b} \alpha \mathbf{c} \mathrm{Cos} \left(\frac{2\pi}{\lambda} (\mathbf{Z} - \mathbf{c}\mathbf{t})\right)} \mathrm{at} \, \mathbf{R} = \mathbf{R}\_1. \tag{14}$$

If (R, Z, U, W) and (r, z, u, w) are, respectively, the coordinates and velocities in the laboratory and wave frame, then transformations from the laboratory frame to wave frame for a steady problem are [48]:

$$\begin{array}{l} \mathbf{r} = \mathbf{R}, \mathbf{z} = \mathbf{Z} - \mathbf{c}t, \mathbf{p}(\mathbf{r}, \mathbf{z}) = \mathbf{P}(\mathbf{R}, \mathbf{Z} - \mathbf{c}t), \mathbf{u}(\mathbf{r}, \mathbf{z}) = \mathbf{U}(\mathbf{R}, \mathbf{Z} - \mathbf{c}t), \\ \mathbf{w}(\mathbf{r}, \mathbf{z}) = \mathbf{W}(\mathbf{R}, \mathbf{Z} - \mathbf{c}t) - \mathbf{c}. \end{array} \tag{15}$$

We introduce the following dimensionless quantities in the wave frame as [29,31]:

r = <sup>r</sup> a2 , z = <sup>z</sup> <sup>λ</sup> , <sup>δ</sup> <sup>=</sup> a2 <sup>λ</sup> ,r1 <sup>=</sup> r1 a2 <sup>=</sup> <sup>ξ</sup>, r2 <sup>=</sup> r2 a2 , u = <sup>λ</sup><sup>u</sup> a2c , <sup>w</sup> <sup>=</sup> <sup>w</sup> <sup>c</sup> , <sup>t</sup> <sup>=</sup> ct λ , φ = <sup>φ</sup> H0a2 ,ψ = <sup>ψ</sup> a2c , Hr <sup>=</sup> <sup>−</sup><sup>δ</sup> r ∂φ <sup>∂</sup><sup>z</sup> , Hz <sup>=</sup> <sup>1</sup> r ∂φ <sup>∂</sup><sup>r</sup> ,u <sup>=</sup> <sup>−</sup><sup>δ</sup> r ∂ψ . ∂ ¯ z , ε = <sup>b</sup> a2 , w = <sup>1</sup> r ∂ψ <sup>∂</sup><sup>r</sup> , <sup>p</sup> <sup>=</sup> a2 <sup>n</sup><sup>+</sup>1p cnλ<sup>k</sup> , <sup>θ</sup> <sup>=</sup> <sup>T</sup>−T1 T0−T1 , <sup>E</sup> = <sup>−</sup> <sup>E</sup> cH0μ<sup>ˆ</sup> . (16)

A non-dimensional governing model for the aforementioned quantities along with long wavelength and creeping Stokesian flow approach is:

$$\begin{split} \frac{\partial \mathbf{p}}{\partial \mathbf{r}} + \frac{1}{\mathbf{r}} \Big( -\frac{1}{\mathbf{r}^2} \frac{\partial \psi}{\partial \mathbf{r}} + \frac{1}{\mathbf{r}} \frac{\partial^2 \psi}{\partial \mathbf{r}^2} \Big)^2 + 2 \Big( -\frac{1}{\mathbf{r}^2} \frac{\partial \psi}{\partial \mathbf{r}} + \frac{1}{\mathbf{r}} \frac{\partial^2 \psi}{\partial \mathbf{r}^2} \Big) \Big( \frac{2}{\mathbf{r}^3} \frac{\partial \psi}{\partial \mathbf{r}} - \frac{2}{\mathbf{r}^2} \frac{\partial^2 \psi}{\partial \mathbf{r}^2} + \frac{1}{\mathbf{r}} \frac{\partial^3 \psi}{\partial \mathbf{r}^3} \Big) - \\ \mathbf{M}^2 \Big( \mathbf{E} - \frac{1}{\mathbf{r}} \frac{\partial \psi}{\partial \mathbf{r}} \Big) - \mathbf{A}\_1 \text{Gr} \theta = 0, \end{split} \tag{17}$$

$$\frac{\partial \mathbf{p}}{\partial \mathbf{r}} = 0.\tag{18}$$

Equation (17) is simplified in the form:

$$-\frac{1}{r^2} \left( -\frac{1}{r^2} \frac{d\psi}{dr} + \frac{1}{r} \frac{d^2 \psi}{dr^2} \right)^2 + \frac{2}{r} \left( -\frac{1}{r^2} \frac{d\psi}{dr} + \frac{1}{r} \frac{d^2 \psi}{dr^2} \right) \left( \frac{2}{r^3} \frac{d\psi}{dr} - \frac{2}{r^2} \frac{d^2 \psi}{dr^2} + \frac{1}{r} \frac{d^3 \psi}{dr^3} \right)$$

$$+ 2 \left( -\frac{1}{r^2} \frac{d\psi}{dr} + \frac{1}{r} \frac{d^2 \psi}{dr^2} \right) \left( -\frac{6}{r^4} \frac{d\psi}{dr} + \frac{6}{r^3} \frac{d^2 \psi}{dr^2} - \frac{3}{r^2} \frac{d^3 \psi}{dr^3} + \frac{1}{r} \frac{d^4 \psi}{dr^4} \right) \tag{19}$$

$$+ 2 \left( \frac{2}{r^3} \frac{d\psi}{dr} - \frac{2}{r^2} \frac{d^2 \psi}{dr^2} + \frac{1}{r} \frac{d^3 \psi}{dr^3} \right)^2 + M^2 \left( -\frac{1}{r^2} \frac{d\psi}{dr} + \frac{1}{r} \frac{d^2 \psi}{dr^2} \right) - A\_1 G r \frac{d\mathcal{O}}{dr} = 0,$$

$$110 \qquad r^2 \Omega \qquad Q.$$

$$\frac{1}{\mathbf{r}}\frac{\partial\theta}{\partial\mathbf{r}} + \frac{\partial^2\theta}{\partial\mathbf{r}^2} + \frac{\Omega}{\mathbf{A}\_2} = 0,\tag{20}$$

$$\mathbf{E} - \frac{1}{\mathbf{r}} \frac{\partial \psi}{\partial \mathbf{r}} - \frac{1}{\mathbf{R}\_{\mathbf{m}}} \left( -\frac{1}{\mathbf{r}^2} \frac{\partial \phi}{\partial \mathbf{r}} + \frac{1}{\mathbf{r}} \frac{\partial^2 \phi}{\partial \mathbf{r}^2} \right) = 0. \tag{21}$$

where, bar notation is ignored.

Corresponding boundary conditions are listed as:

$$\begin{aligned} \psi(\mathbf{r}) &= -\frac{\mathbb{F}}{2}, \; \frac{1}{\mathbf{r}} \frac{\partial \psi}{\partial \mathbf{r}} = -1, \; \phi(\mathbf{r}) = 0, \; \theta(\mathbf{r}) = 1, \; \text{at } \mathbf{r} = \mathbf{r}\_1 = \boldsymbol{\xi}, \\\ \psi(\mathbf{r}) &= \frac{\mathbb{F}}{2}, \; \frac{1}{\mathbf{r}} \frac{\partial \psi}{\partial \mathbf{r}} = -1 - 2\pi\varepsilon a\delta \cos(2\pi\mathbf{z}), \; \phi(\mathbf{r}) = 1, \\\ \theta(\mathbf{r}) &= 0, \; \text{at } \mathbf{r} = \mathbf{r}\_2 = 1 + \varepsilon \cos(2\pi\mathbf{z}). \end{aligned} \tag{22}$$

In the above expressions, *u* and *w* denotes *r*- and *z*-components of velocity within the wave frame, respectively. Emerging parameters in the above model are expressed as [48,49]:

$$\begin{split} \text{M}^{2} &= \text{Re}\mathbf{S}^{2} \mathbf{R}\_{\text{m}}, \text{Re} = \frac{\mathbf{a}\_{2}^{\text{n}} \rho\_{\text{t}}}{\text{k} \mathbf{c}^{n-2}}, \mathbf{R}\_{\text{m}} = \sigma \text{fl} \mathbf{a}\_{2} \mathbf{c}, \mathbf{S} = \frac{\text{H}\_{0}}{\text{c}} \sqrt{\frac{\mu}{\rho\_{\text{t}}}}, \\ \text{p}\_{\text{m}} &= \text{p} + \frac{1}{2} \text{Re}\delta \frac{\rho (\text{H})^{2}}{\text{c}^{2} \rho\_{\text{t}}}, \text{ } \mathcal{Q} = \frac{\text{Q} \mathbf{a}\_{2}^{2}}{\text{x}\_{1} (\text{T}\_{0} - \text{T}\_{1})}, \text{ Gr} = \frac{(\rho \emptyset)\_{\text{t}} (\text{T}\_{0} - \text{T}\_{1}) \mathbf{a}\_{2}^{n+1}}{\text{k} \mathbf{c}^{n}}. \end{split} \tag{23}$$

where as

$$\begin{split} \mathcal{A}\_{1} &= \left(1 - \phi\_{1} - \phi\_{2}\right) + \phi\_{1} \big(\frac{(\rho\beta)\_{\mathbf{v}\_{1}}}{(\rho\beta)\_{\mathbf{f}}}\big) + \phi\_{2} \frac{(\rho\beta)\_{\mathbf{v}\_{2}}}{(\rho\beta)\_{\mathbf{f}}},\\ \mathcal{A}\_{2} &= \frac{\kappa\_{\mathbf{s}\_{2}} + (s-1)\kappa\_{\mathbf{h}\mathbf{f}} - (s-1)\phi\_{2}(\kappa\_{\mathbf{h}\mathbf{f}} - \kappa\_{\mathbf{s}\_{2}})}{\kappa\_{\mathbf{s}\_{2}} + (s-1)\kappa\_{\mathbf{h}\mathbf{f}} + \phi\_{2}(\kappa\_{\mathbf{h}\mathbf{f}} - \kappa\_{\mathbf{s}\_{2}})} \frac{\kappa\_{\mathbf{s}\_{1}} + (s-1)\kappa\_{\mathbf{f}} - (s-1)\phi\_{1}(\kappa\_{\mathbf{f}} - \kappa\_{\mathbf{s}\_{1}})}{\kappa\_{\mathbf{s}\_{1}} + (s-1)\kappa\_{\mathbf{f}} + \phi\_{1}(\kappa\_{\mathbf{f}} - \kappa\_{\mathbf{s}\_{1}})}. \end{split} \tag{24}$$

Moreover, the pressure gradient can be achieved from the following relation:

$$F = \int\_{r\_1}^{r\_2} rw dr = \int\_{r\_1}^{r\_2} \frac{\partial \psi}{\partial r} dr \tag{25}$$

where, *F* is the volumetric rate of flow inside the wave frame. Now, non-dimensional mean flow rate *Q* into the laboratory frame assuming the transformations of Equations (16) is:

$$Q = F + \frac{1}{2} \left( 1 - \xi^2 + \frac{\varepsilon^2}{2} \right). \tag{26}$$

Pressure rise per wavelength is calculated utilizing Equation (25) as:

$$
\Delta P = \int\_0^1 \frac{dp}{dz} dz.\tag{27}
$$

All variables and parameter are defined in Appendix A.

#### **3. Methodology and Convergence of HAM Solutions**

### *3.1. Methodology*

The dimensionless governing model containing Equations (19)–(21) under the associated boundary conditions (22) is analyzed by employing homotopy analysis method. For this, the initial guesses are ψ0(*r*), φ0(*r*) and θ0(*r*) and linear operators are chosen in the subsequent manner as:

$$L\_1(\psi) = \psi^{(iv)},\ L\_2(\phi) = \phi^{\prime\prime},\\ L\_3(\theta) = \theta^{\prime}.\tag{28}$$

And

$$L\_1 \left(\mathbb{C}\_1 + \mathbb{C}\_2 r + \mathbb{C}\_3 r^2 + \mathbb{C}\_4 r^3\right) = L\_2 \left(\mathbb{C}\_5 + \mathbb{C}\_6 r\right) = L\_3 \left(\mathbb{C}\_7 + \mathbb{C}\_8 r\right) = 0,\tag{29}$$

where, *C1-C8* represents constants while *h*1, *h*<sup>2</sup> and *h*<sup>3</sup> being auxiliary parameter which plays a key role in the frame of HAM, since the convergence of solutions strongly depends on *h*. Now, for embedding parameter γ ∈ [0, 1] and non-zero auxiliary parameters, the problem under study can be constructed in the following manner [29,30]:

Zeroth-order deformation problem:

$$
\hat{\eta}\_1(1-\gamma)L\_1[\psi(r,\gamma)-\psi\_0(r)]=\gamma h\_1 N\_1[\psi(r,\gamma),\phi(r,\gamma),\theta(r,\gamma)],\tag{30}
$$

$$\frac{1}{2}(1-\gamma)L\_2[\phi(r,\gamma)-\phi\_0(r)]=\gamma h\_2 N\_2[\psi(r,\gamma),\phi(r,\gamma),\theta(r,\gamma)],\tag{31}$$

$$
\hat{\rho}\_1(1-\gamma)L\_3[\theta(r,\gamma)-\theta\_0(r)]=\gamma h\_3 N\_3[\psi(r,\gamma),\phi(r,\gamma),\theta(r,\gamma)].\tag{32}
$$

and,

$$\begin{aligned} \text{At } r &= r\_1 = \xi, \psi(r\_1, \gamma) = -\frac{\mathbb{F}}{2}, \frac{1}{r\_1} \psi'(r\_1, \gamma) = -1, \phi(r\_1, \gamma) = 0, \theta(r\_1, \gamma) = 1. \\ \text{At } r &= r\_2 = 1 + \varepsilon \cos(2\pi z), \psi(r\_2, \gamma) = \frac{\mathbb{F}}{2}, \frac{1}{r\_2} \psi'(r\_2, \gamma) = -1 - 2\pi\varepsilon a\delta \cos(2\pi z), \\ \phi(r\_2, \gamma) &= 1, \theta(r\_2, \gamma) = 0. \end{aligned} \tag{33}$$

On the basis of selected linear operator, auxiliary parameter and initial guesses, the m th order solution series is constructed as:

$$L\_1[\psi\_m(r,\gamma) - \chi\_m \psi\_{m-1}(r,\gamma)] = h\_1 R\_m^1(r,\gamma),\tag{34}$$

$$L\_2[\phi\_m(r,\gamma) - \chi\_m \phi\_{m-1}(r,\gamma)] = h\_2 R\_m^2(r,\gamma),\tag{35}$$

$$L\_3[\theta\_m(r,\gamma) - \chi\_m \theta\_{m-1}(r,\gamma)] = h\_3 R\_m^3(r,\gamma). \tag{36}$$

The boundary conditions are:

$$\text{at } r = r\_1 = \xi, \psi\_m(r\_1, \gamma) = 0, \frac{1}{r\_1} \psi'\_m(r\_1, \gamma) = 0, \phi\_m(r\_1, \gamma) = 0, \theta\_m(r\_1, \gamma) = 0.$$

$$\text{at } r = r\_2 = 1 + \varepsilon \cos(2\pi z), \psi\_m(r\_2, \gamma) = 0, \frac{1}{r\_2} \psi'\_m(r\_2, \gamma) = 0, \phi\_m(r\_2, \gamma) = 0,\tag{37}$$

$$\theta\_m(r\_2, \gamma) = 0.$$

where the auxiliary parameter is found by plotting *h*-curves while χ*<sup>m</sup>* is defined as:

$$\chi\_m = \begin{array}{c} 0, & m \le 1 \\ 1, & m > 1 \end{array} \Big| 1$$

Therefore, we can write:

$$\psi\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) = \mathbb{V}\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) + \mathbb{C}\_{1} + \mathbb{C}\_{2}r + \mathbb{C}\_{3}r^{2} + \mathbb{C}\_{4}r^{3},$$

$$\phi\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) = \Phi\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) + \mathbb{C}\mathfrak{s} + \mathbb{C}\_{6}r,\tag{38}$$

$$\theta\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) = \Theta\_{\mathfrak{m}}(r,\boldsymbol{\gamma}) + \mathbb{C}\boldsymbol{\gamma} + \mathbb{C}\mathfrak{s}r.$$

The particular solutions Ψ*m*(*r*, γ), Φ*m*(*r*, γ) and Θ*m*(*r*, γ) are obtained using a symbolic software mathematica while the constants are determined from the defined boundary conditions.
