MHD flow and heat transfer of couple stress fluid over an oscillatory stretching sheet with heat source/sink in porous medium

==Abstract==

This article presents the MHD flow and heat transfer of a couple stress fluid over oscillatory stretching sheet embedded in a porous medium in the presence of heat source/sink. The unsteady flow problem is reduced to two coupled partial differential equations using dimensionless variables. Homotopy analysis method is employed to obtain the solution of these equations. Based on the solution of these equations an extensive analysis is performed to investigate the effects of various flow parameters on velocity and temperature fields, skin friction coefficient and Nusselt number. It is found that the presence of couple stress in viscous fluid increases the amplitude of oscillations in velocity and skin friction coefficient. It is also noticed that temperature increases by increasing heat source parameter. Moreover, the numerical values of local Nusselt number are calculated and shown in tabular form. It is found that the Nusselt number increases by increasing Prandtl number while it decreases by increasing couple stress parameter.

==Keywords==

Boundary layer flow; Couple stress fluid; Oscillatory stretching sheet; Homotopy analysis method

==1. Introduction==

It is well established fact that Navier–Stokes equations are not adequate choice in describing the flows of non-Newtonian fluids. Non-Newtonian fluids frequently occur in many industrial and technological processes. If the flow of non-Newtonian fluid is considered over a stretching sheet then it become relevant to many industrial applications such a extrusion of plastic sheets, paper production, crystal growing and glass blowing, drawing of copper wires, continuous stretching of plastic films. The pioneering work on viscous flow due to continuous moving sheet was performed by Sakiadis [[#b0005|[1]]]. His analysis was extended for stretching sheet by Crane [[#b0010|[2]]]. Crane’s problem was extended by Rajagopal et al. [[#b0015|[3]]] for viscoelastic fluid. Gupta and Gupta [[#b0020|[4]]] extended the analysis of Crane by considering suction/blowing at the sheet surface. Dandapat and Gupta [[#b0025|[5]]] discussed the heat transfer characteristics of viscoelastic fluid over a stretching sheet. A study of uniqueness of boundary value problem arising due to a stretching boundary in a viscous fluid was performed by McLeod and Rajagopal [[#b0030|[6]]]. Rollins and Vajravelu [[#b0035|[7]]] discussed heat transfer characteristics of second-order fluid over a continuos stretching surface with internal heat generation or absorption. The stretching sheet problem for viscoelastic Oldroyd-B model was extended by Rajagopal and Bhatnagar [[#b0040|[8]]]. Cortell [[#b0045|[9]]] presented a study of flow and heat transfer of hydromagnetic second-order fluid past a stretching sheet. In another paper [[#b0050|[10]]], Cortell included magnetic and suction effects while discussing the flow and heat transfer of second grade fluid due to stretching sheet. The study of stagnation-point flow of micropolar fluid toward a stretching sheet was initiated by Nazar et al. [[#b0055|[11]]]. Combined forced and natural convection flow of micropolar fluid over a stretching sheet in the vicinity of stagnation point was discussed by Ishak et al. [[#b0060|[12]]]. Ariel [[#b0065|[13]]] put forward a study of flow of second grade fluid past a radially stretching sheet. Ariel work was extended by Hayat and Sajid [[#b0070|[14]]] to include the heat transfer characteristics. Some other contribution by Hayat and his co-workers regarding heat transfer analysis over a stretching surface can be found in references [[#b0075|[15]]], [[#b0080|[16]]], [[#b0085|[17]]] and [[#b0090|[18]]]. The investigation of non-similar flow of third order fluid over a stretching sheet was due to Sajid and Hayat [[#b0095|[19]]]. Hayat et al. [[#b0100|[20]]] analyzed combined heat and mass transfer characteristics in boundary layer flow of viscoelastic fluid over linearly stretching vertical surface by using homotopy analysis method. Later on, Turkyilmazoglu [[#b0105|[21]]] computed exact solution of hydromagnetic flow of viscoelastic fluid in the presence of Soret and Dufour effects. Turkyilmazoglu [[#b0110|[22]]] also developed closed-form solutions for three dimensional flow and heat transfer of a electrically conducting viscoelastic fluid over a porous stretching/shrinking sheet. The multiplicity of the solutions in boundary layer flows over stretching/shrinking sheet is also addressed by Turkyilmazoglu. In fact, it was shown by him that multiple solutions exist for two-dimensional flow over stretching sheet and in the case of shrinking sheet triple solutions can be detected [[#b0115|[23]]], [[#b0120|[24]]] and [[#b0125|[25]]]. In another paper, Turkyilmazoglu [[#b0130|[26]]] computed exponential series type solution in the hydromagnetic and heat transfer in two–three dimensional flow over stretching or shrinking surface. According to him exponential series type solution exists with finite and infinite terms in the expansion depending upon prevailing parameters in the governing equations.

To best of our knowledge, Wang [[#b0135|[27]]] initiated the study of viscous fluid due to oscillatory stretching sheet. His work was extended by Saddipa et al. [[#b0140|[28]]] by using constitutive equation of Walter’s B model. Later on, Abbas et al. [[#b0145|[29]]] complemented the Wang ’s analysis by including the effects of slip and heat transfer. They also utilized the constitutive relation of second order fluid [[#b0150|[30]]] to study Wang’s problem. Zheng et al. [[#b0155|[31]]] provided an analysis of unsteady heat and mass transfer in MHD flow over an oscillatory stretching sheet with Soret and Dufour effects. They have employed homotopy analysis method and delineated the effects of various emerging parameters on flow and heat transfer characteristics. Their study indicates a decrease in velocity of fluid by increasing unsteady parameter, magnetic parameter and suction/injection parameter. They also found an increase in thermal boundary layer thickness with magnetic parameter and unsteady parameter. Ali et al. [[#b0160|[32]]] analyzed the unsteady flow of Jeffrey fluid over an oscillatory stretching surface. Khan et al. [[#b0165|[33]]] used an implicit finite difference scheme to analyze the heat transfer effects in the hydromagnetic flow of viscoelastic fluid (second grade fluid) over an oscillatory stretching surface. However, so far no one has considered the study of unsteady boundary layer flow of a couple stress fluid through a porous medium over an oscillatory stretching surface in the presence of heat generation or absorption. Hence we feel appropriate to consider the study of such flows, as such type of flows may find a number of industrial applications. The theory of couple stress fluid was developed by Stokes [[#b0170|[34]]] as a simple generalization of classical viscous theory that sustains couple stresses and body couples. The main effects of couple stress fluid are to introduce a size dependent effect which is not present in classical viscous theories. Blood, lubricants containing small amount of additives, electro-rheological and synthetic fluids are examples of couple stress fluid [[#b0175|[35]]]. Couple stress theory has been successfully applied to many problems in biomechanics and lubrications area by Srivastava [[#b0180|[36]]], El-Shehawey and Mekheimer [[#b0185|[37]]], Pal et al. [[#b0190|[38]]], Chiang et al. [[#b0195|[39]]], Naduvinamani et al. [[#b0200|[40]]], Jian and Chen et al. [[#b0205|[41]]] and Lu and Lin [[#b0210|[42]]].

The presentation of the paper is as follows. Flow analysis is presented in Section [[#s0010|2]]. Homotopy analysis method is employed for solution of governing equations in Section [[#s0015|3]]. Convergence of solution is discussed in Section [[#s0020|4]]. Quantitative analysis of results is presented in Section [[#s0025|5]]. Final remarks are presented at the end of the paper in Section [[#s0030|6]].

<span id='s0010'></span>
==2. Flow analysis==

Let us consider the unsteady two-dimensional magnetohydrodynamics flow of an incompressible fluid over an oscillatory stretching sheet. In a rectangular coordinate system <math display="inline">\left(\overline{x}\mbox{,}\overline{y}\right)</math>, the sheet is assumed to coincide with <math display="inline">\overline{y}=0</math> and flow takes place in the semi-infinite porous space <math display="inline">\left(\overline{y}>0\right)</math> of constant permeability (see [[#f0005|Fig. 1]]). A constant magnetic field of strength <math display="inline">B_0</math> is applied perpendicular to the stretching surface. It is assumed that the sheet is stretched with velocity <math display="inline">u_{\omega }=b\overline{x}sin\omega t</math> in which ''b  '' is the stretching rate and <math display="inline">\omega </math> denotes the angular velocity. Moreover, sheet is maintained at uniform temperature <math display="inline">T_w</math><math display="inline">\left(>T_{\infty }\right)</math>, where <math display="inline">T_{\infty }</math> is the temperature of ambient fluid. Now problem is to determine the velocity and temperature fields inside the fluid satisfying the boundary conditions at the wall and far away from the wall. To this end, we apply boundary layer approximation on continuity and momentum equations and write

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\frac{\partial u}{\partial \overline{x}}+\frac{\partial v}{\partial \overline{y}}=</math><math>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|}

<span id='e0010'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial \overline{x}}+</math><math>v\frac{\partial u}{\partial \overline{y}}=\nu \frac{{\partial }^2u}{\partial {\overline{y}}^2}-</math><math>\frac{{\eta }_0}{\rho }\left[\frac{{\partial }^4u}{\partial {\overline{y}}^4}\right]-</math><math>\frac{\sigma B_0^2}{\rho }u-\frac{\upsilon \phi }{k}u\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|}

<span id='e0015'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\frac{\partial T}{\partial t}+u\frac{\partial T}{\partial \overline{x}}+</math><math>v\frac{\partial T}{\partial \overline{y}}=\frac{k}{\rho c_p}\frac{{\partial }^2T}{\partial {\overline{y}}^2}+</math><math>\frac{Q}{\rho c_p}\left(T-T_{\infty }\right)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|}

where ''u'' and ''v  '' are velocity component along <math display="inline">\overline{x}</math> and <math display="inline">\overline{y}</math>-directions, respectively, <math display="inline">\nu </math> is the kinematic viscosity, <math display="inline">\rho </math> is the density, <math display="inline">{\eta }_0</math> is the material constant for the couple stress fluid, ''Q  '' is the heat source coefficient, <math display="inline">\sigma </math> is electric conductivity and <math display="inline">B_0</math> is the strength of constant applied magnetic field. The flow is subjected to the following boundary conditions

<span id='e0020'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>u=u_{\omega }=b\overline{x}sin\omega t\mbox{,}\quad v=</math><math>0\mbox{,}\quad \frac{{\partial }^2u}{\partial {\overline{y}}^2}=</math><math>0\mbox{,}\quad T=T_w\quad at\quad \overline{y}=0\mbox{,}\quad t>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|}

<span id='e0025'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>u=0\mbox{,}\quad \frac{\partial u}{\partial \overline{y}}=</math><math>0\mbox{,}\quad T\rightarrow T_{\infty }\quad \quad \quad at\quad \quad \quad \overline{y}\rightarrow \infty \mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|}

The first two boundary conditions in [[#e0020|(4)]] result from the assumption of no-slip at the surface while the third condition is consequence of the assumption of vanishing of couple stresses at the wall. Since flow is in unbounded domain that is why augmented boundary condition is imposed at <math display="inline">y\rightarrow \infty </math>. To non-dimensionalize the flow problem, we use the dimensionless variables  [[#b0150|[30]]] and [[#b0155|[31]]]

<span id='e0030'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>y=\sqrt{\frac{b}{\nu }}\overline{y}\mbox{,}\quad \tau =</math><math>t\omega \mbox{,}\quad u=b\overline{x}f_y\left(y\mbox{,}\tau \right)\mbox{,}\quad v=</math><math>-\sqrt{\nu b}f\left(y\mbox{,}\tau \right)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
|}

<span id='e0035'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\theta (y\mbox{,}\tau )=\frac{T-T_{\infty }}{T_w-T_{\infty }}\mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
|}

In view of above transformations, the continuity equation is identically satisfied and Eqs.  [[#e0010|(2)]] and [[#e0015|(3)]] are transformed as follows

<span id='e0040'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_{yyy}-{Sf}_{y\tau }-f_y^2+{ff}_{yy}-\beta f_y-{Kf}_{yyyyy}=</math><math>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (8)
|}

<span id='e0045'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\theta }_{yy}+Pr\left(f{\theta }_y-S{\theta }_{\tau }+\right. </math><math>\left. \lambda \theta \right)=0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
|}

where the subscript denotes differentiation with respect to ''y''. The boundary conditions  [[#e0020|(4)]] and [[#e0025|(5)]] take the following form

<span id='e0050'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_y\left(0\mbox{,}\tau \right)=sin\tau \mbox{,}\quad f\left(0\mbox{,}\tau \right)=</math><math>0\mbox{,}\quad f_{yyy}\left(0\mbox{,}\tau \right)=</math><math>0\mbox{,}\quad \theta (0\mbox{,}\tau )=1\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
|}

<span id='e0055'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_y\left(\infty \mbox{,}\tau \right)=0\mbox{,}\quad f_{yy}\left(\infty \mbox{,}\tau \right)=</math><math>0\mbox{,}\quad \theta (\infty \mbox{,}\tau )=0\mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
|}

In above equations <math display="inline">K={\eta }_0b/\rho {\upsilon }^2</math> is couple stress parameter, <math display="inline">S=\omega /b</math> is unsteady parameter, <math display="inline">\beta =\sigma B_0^2/\rho b+\upsilon \phi /kb</math> is the magneto-porous parameter, <math display="inline">Pr=\mu c_p/k</math> is the Prandtl number and <math display="inline">\lambda =Q/c\rho c_p</math> is the heat source <math display="inline">\left(\lambda >0\right)</math> or sink <math display="inline">\left(\lambda <0\right)</math> parameter. The parameter <math display="inline">K={\eta }_0\rho b/{\mu }^2=l^2b\rho /\mu </math> is called couple stress parameter. It has been pointed out by Stokes [[#b0170|[34]]] that the effects of couple stress are quite large for large values of <math display="inline">K=l^2b\rho /\mu </math>, where <math display="inline">b\rho /\mu </math> is a typical length scale associated with the flow and <math display="inline">l=\sqrt{\eta /\mu }</math> is a material coefficient which is a function of molecular dimensions of the liquid. The parameter ''l'' will be greatly varied for different liquids. For instant, the length of polymer chain may be a million time the diameter of a water molecule. This is the size dependent effect associated with the couple stress fluid theory. For small value of ''K'' the effects of couple stress diminish and fluid behaves like a Newtonian fluid.

<span id='f0005'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr1.jpg|center|311px|Geometry of the problem.]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 1.

Geometry of the problem.

</span>
|}

The quantities of interest are skin friction coefficient and local Nusselt number which are defined as

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>C_f=\frac{{\tau }_w}{\rho u_w^2}\mbox{,}\quad {Nu}_x=</math><math>\frac{\overline{x}q_w}{k(T_w-T_{\infty })}\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (12)
|}

where

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\tau }_w=\mu {\left(\frac{\partial u}{\partial \overline{y}}\right)}_{\overline{y}=0}-</math><math>\frac{{\eta }_0}{\rho }{\left(\frac{{\partial }^3u}{\partial {\overline{y}}^3}\right)}_{\overline{y}=0}\mbox{,}\quad \quad \quad \quad \quad q_w=</math><math>-k{\left(\frac{\partial T}{\partial \overline{y}}\right)}_{\overline{y}=0}\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
|}

are skin friction and heat transfer at the wall, respectively. In view of transformations [[#e0030|(6)]] and [[#e0035|(7)]], we get

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{Re}_x^{1/2}C_f=f_{yy}(0\mbox{,}\tau )-{Kf}_{yyyy}(0\mbox{,}\tau )\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (14)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{Re}_x^{-1/2}{Nu}_x=-{\theta }_y(0\mbox{,}\tau )\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (15)
|}

where <math display="inline">{Re}_x=u_w\overline{x}/\nu </math> is the local Reynolds number.

<span id='s0015'></span>
==3. Solution by homotopy analysis method==

In many physical problems, the resulted differential equations are highly nonlinear nature. The computations of analytical or numerical solutions of such equations have set a challenge for the researchers. Homotopy analysis method (HAM) is one the best analytical methods to compute the series solution of such nonlinear partial as well as ordinary differential equations. This method tackles for strongly nonlinear problems even if a given nonlinear problem does not contain any small/large parameter. This method gives great freedom to choose and adjust the convergence region and rate of approximation. The homotopy analysis method is advantageous over routine numerical methods in the sense that it is free of rounding off errors which arise due to discretization procedure. Moreover, it does not require large computer memory and time [[#b0215|[43]]]. This method was successfully applied by researchers to solve various nonlinear problems in various disciplines of science and engineering [[#b0220|[44]]], [[#b0225|[45]]], [[#b0230|[46]]], [[#b0235|[47]]], [[#b0240|[48]]] and [[#b0245|[49]]]. It has been emphasized by Liao that the series solution obtained by HAM is purely analytic and he demonstrated this fact in number of articles for instance [[#b0250|[50]]], [[#b0255|[51]]] and [[#b0260|[52]]]. Moreover, it has been proved by Turkyilmazoglu [[#b0265|[53]]] that HAM can handle any type of nonlinearity by a suitable HAM formulation. For the series solution of Eqs. [[#e0040|(8)]] and [[#e0045|(9)]] subject to the boundary conditions [[#e0050|(10)]] and [[#e0055|(11)]], we employ homotopy analysis method (HAM). According to HAM the velocity and temperature fields can be represented as

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f(y\mbox{,}\tau )=a_{0\mbox{,}0}^0+\sum_{n=0}^{\infty }\sum_{k=0}^{\infty }\sum_{j=0}^{\infty }a_{n\mbox{,}k}^jy^ksin(j\tau )exp(-</math><math>ny)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (16)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\theta (y\mbox{,}\tau )=b_{0\mbox{,}0}^0+\sum_{n=0}^{\infty }\sum_{k=0}^{\infty }\sum_{j=0}^{\infty }b_{n\mbox{,}k}^jy^ksin(j\tau )exp(-</math><math>ny)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (17)
|}

where <math display="inline">a_{n\mbox{,}k}^j</math> and <math display="inline">b_{n\mbox{,}k}^j</math> are the coefficients to be determined. The appropriate initial guesses for <math display="inline">f(y\mbox{,}\tau )</math> and <math display="inline">\theta (y)</math> are

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_0(y\mbox{,}\tau )=\frac{1}{2}sin\tau (3-3exp(-y)-</math><math>yexp(-y))\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (18)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\theta }_0(y)=exp(-y)\mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (19)
|}

Further we choose

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_f(f)=\frac{{\partial }^3f}{\partial y^3}-</math><math>\frac{\partial f}{\partial y}\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | 
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_{\theta }(f)=\frac{{\partial }^2f}{\partial y^2}-</math><math>f\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | 
|}

as auxiliary linear operators having following properties

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_f[C_1+C_2exp(-y)+C_3exp(y)]=0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (20)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_{\theta }[C_4exp(y)+C_5exp(-y)]=0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (21)
|}

where <math display="inline">C_i</math><math display="inline">\left(i=1-5\right)</math> are arbitrary constants. If <math display="inline">h_f</math> and <math display="inline">h_{\theta }</math> denotes the nonzero auxiliary parameters then zeroth order deformation problems are constructed as follows:

<span id='e0110'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>(1-p){\pounds }_f[\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)-</math><math>f_0(y\mbox{,}\tau )={ph}_fN_f[\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)]\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (22)
|}

<span id='e0115'></span>
{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>(1-p){\pounds }_f[\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)-</math><math>{\theta }_0(y\mbox{,}\tau )={ph}_{\theta }N_{\theta }[\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)\mbox{,}\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)]\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (23)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{f}(0\mbox{,}\tau \mbox{;}p)=0\mbox{,}\frac{\partial \overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y}=</math><math>sin\tau \mbox{,}\quad \frac{\partial \overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y}=</math><math>0\mbox{,}\quad \frac{{\partial }^2\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y^2}=</math><math>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (24)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{\theta }(0\mbox{,}\tau \mbox{;}p)=</math><math>1\mbox{,}\quad \overset{\mbox{ˆ}}{\theta }(\infty \mbox{,}\tau \mbox{;}p)=</math><math>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | 
|}

where <math display="inline">p\in [0\mbox{,}1]</math> is an embedding parameter and nonlinear operators <math display="inline">N_f</math> and <math display="inline">N_{\theta }</math> are

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>N_f[\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)]=</math><math>-K\frac{{\partial }^5\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y^5}+</math><math>\frac{{\partial }^3\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y^3}-</math><math>S\frac{{\partial }^2\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y\partial \tau }+</math><math>\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)\frac{{\partial }^2\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y^2}-</math><math>{\left(\frac{\partial \overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y}\right)}^2-</math><math>\beta \frac{\partial \overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial y}\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (25)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>N_{\theta }[\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)\mbox{,}\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)]=</math><math>\frac{{\partial }^2\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)}{\partial y^2}+</math><math>Pr\left(\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)\frac{\partial \overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)}{\partial y}-\right. </math><math>\left. S\frac{\partial \overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)}{\partial \tau }+\right. </math><math>\left. \lambda \overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)\right)\mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (26)
|}

The zeroth-order deformation problems defined above have the following solutions corresponding to <math display="inline">p=0</math> and <math display="inline">p=1</math>

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}0)=f_0(y\mbox{,}\tau )\mbox{,}\quad \overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}1)=</math><math>f(y\mbox{,}\tau )\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (27)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}0)=</math><math>{\theta }_0(y\mbox{,}\tau )\mbox{,}\quad \overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}1)=</math><math>\theta (y\mbox{,}\tau )\mbox{.}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (28)
|}

Expanding <math display="inline">\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)</math> and <math display="inline">\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)</math> in a Taylor’s series with respect to ''p'', we get

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)={\overset{\mbox{ˆ}}{f}}_0(y\mbox{,}\tau )+</math><math>\sum_{m=1}^{\infty }f_m(y\mbox{,}\tau )p^m\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (29)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)=</math><math>{\overset{\mbox{ˆ}}{\theta }}_0(y\mbox{,}\tau )+\sum_{m=1}^{\infty }{\theta }_m(y\mbox{,}\tau )p^m\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (30)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_m(y\mbox{,}\tau )=\frac{1}{m!}\frac{{\partial }^m\overset{\mbox{ˆ}}{f}(y\mbox{,}\tau \mbox{;}p)}{\partial p^m}\mbox{,}\quad {\theta }_m(y\mbox{,}\tau )=</math><math>\frac{1}{m!}\frac{{\partial }^m\overset{\mbox{ˆ}}{\theta }(y\mbox{,}\tau \mbox{;}p)}{\partial p^m}\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (31)
|}

Differentiating zeroth order deformations equations  [[#e0110|(22)]] and [[#e0115|(23)]]''m''-times with respect to ''p  '', then setting <math display="inline">p=0</math> and finally dividing by <math display="inline">m!</math>, we get the following ''m''th-order deformation problem

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_f[{\overset{\mbox{ˆ}}{f}}_m(y\mbox{,}\tau )-</math><math>{\varkappa }_mf_{m-1}(y\mbox{,}\tau )=h_fR_m^f(y\mbox{,}\tau )\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (32)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\pounds }_f[{\overset{\mbox{ˆ}}{\theta }}_m(y\mbox{,}\tau )-</math><math>{\varkappa }_m{\theta }_{m-1}(y\mbox{,}\tau )=h_{\theta }R_m^{\theta }(y\mbox{,}\tau )\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (33)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>f_m(0\mbox{,}\tau \mbox{;}p)=0\mbox{,}{\frac{\partial f_m(y\mbox{,}\tau \mbox{;}0)}{\partial y}}_{y=0}=</math><math>{\frac{\partial f_m(y\mbox{,}\tau \mbox{;}0)}{\partial y}}_{y=\infty }=</math><math>{\frac{\partial f_m^2(y\mbox{,}\tau \mbox{;}0)}{\partial y^2}}_{y=\infty }=</math><math>0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (34)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\theta }_m(0)=0\mbox{,}{\theta }_m(\infty )=0\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (35)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>R_m^f(y\mbox{,}\tau )=-K\frac{{\partial }^5f_{m-1}}{\partial y^5}+</math><math>\frac{{\partial }^3f_{m-1}}{\partial y^3}-S\frac{{\partial }^2f_{m-1}}{\partial y\partial \tau }-</math><math>\beta \frac{\partial f_{m-1}}{\partial y}+\sum_{k=0}^{m-1}\left(f_{m-1-k}\frac{{\partial }^2f_k}{\partial y^2}-\right. </math><math>\left. \frac{f_{m-k-1}}{\partial y}\frac{\partial f_k}{\partial y}\right)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (36)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>R_m^{\theta }(y\mbox{,}\tau )=\frac{{\partial }^2{\theta }_{m-1}}{\partial y^2}-</math><math>Pr\left[S\frac{\partial {\theta }_{m-1}}{\partial \tau }\right]+</math><math>Pr\left(\sum_{k=0}^{m-1}\left[f_k\frac{\partial {\theta }_{m-1-k}}{\partial y}\right]+\right. </math><math>\left. \lambda {\theta }_{m-1}\right)\mbox{,}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (37)
|}

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\varkappa }_m=\begin{array}{ll}
0\mbox{,} & m\leqslant 1\mbox{,}\\
1\mbox{,} & m>1\mbox{,}
\end{array}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (38)
|}

<span id='s0020'></span>
==4. Convergence of HAM solution==

It is well established fact that the convergence of the HAM solution is largely dependent on the proper choice of the auxiliary parameters <math display="inline">h_f</math> and <math display="inline">h_{\theta }</math>. For a particular set of parameters the convergence region can be obtained by plotting ''h''-curves.  [[#f0010|Figure 2]] and [[#f0015|Figure 3]] present two such curves showing the plausible values of <math display="inline">h_f</math> and <math display="inline">h_{\theta }</math> for a given set of parameters. We note from these figures that for convergent solution <math display="inline">-1.2\leqslant h_f<-0.2</math> and <math display="inline">-1.5\leqslant h_{\theta }<-0.2</math>. The purpose of showing such curves is just to emphasize that in principle for any physical choice of parameters of the problem, a convergent solution can be obtained. The values of <math display="inline">f^{{''}}(0\mbox{,}\tau )</math> obtained by HAM in the present paper are also compared with the corresponding values reported in Refs.  [[#b0100|[20]]], [[#b0105|[21]]], [[#b0150|[30]]] and [[#b0155|[31]]] in  [[#t0005|Table 1]] and [[#t0010|Table 2]]. Both tables show that our results are in excellent agreement with the existing results in the literature.

<span id='f0010'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr2.jpg|center|311px|The h-curve for velocity at 7th order of approximation: with K=0.1,β=0.2,S=0.1 ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 2.

The ''h  ''-curve for velocity at 7th order of approximation: with <math display="inline">K=0.1\mbox{,}\beta =0.2\mbox{,}S=0.1</math> and <math display="inline">\tau =0.5\pi </math>.

</span>
|}

<span id='f0015'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr3.jpg|center|309px|The h-curve for temperature at 7th order of approximation with ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 3.

The ''h  ''-curve for temperature at 7th order of approximation with <math display="inline">K=0.1\mbox{,}S=0.1\mbox{,}\beta =0.2\mbox{,}\lambda =</math><math>0.5\mbox{,}Pr=0.1</math> and <math display="inline">\tau =0.5\pi </math>.

</span>
|}

<span id='t0005'></span>

{| class="wikitable" style="min-width: 60%;margin-left: auto; margin-right: auto;"
|+

Table 1.

Comparison of values of <math display="inline">f^{{''}}(0\mbox{,}\tau )</math> with Refs. [[#b0150|[30]]] with <math display="inline">K=0</math> and [[#b0155|[31]]].

|-

! ''S''
! <math display="inline">\beta </math>
! <math display="inline">\tau </math>
! Zheng et al. [[#b0155|[31]]]
! Abbas et al. [[#b0150|[30]]]
! Present results
|-

| <math display="inline">1.0</math>
| 12
| <math display="inline">1.5\pi </math>
| <math display="inline">11.678565</math>
| <math display="inline">11.678656</math>
| <math display="inline">11.678565</math>
|-

| 
| 
| <math display="inline">5.5\pi </math>
| <math display="inline">11.678706</math>
| <math display="inline">11.678707</math>
| <math display="inline">11.678706</math>
|-

| 
| 
| <math display="inline">9.5\pi </math>
| <math display="inline">11.678656</math>
| <math display="inline">11.678656</math>
| <math display="inline">11.678656</math>
|}

<span id='t0010'></span>

{| class="wikitable" style="min-width: 60%;margin-left: auto; margin-right: auto;"
|+

Table 2.

Comparison of values of <math display="inline">f^{{''}}(0\mbox{,}\tau )</math> for different values of <math display="inline">\beta </math> when <math display="inline">K=0\mbox{,}S=0</math> and <math display="inline">\tau =\pi /2</math>.

|-

! <math display="inline">\beta </math>
! Hayat et al. [[#b0100|[20]]]
! Turkyilmazoglu [[#b0105|[21]]]
! Present results
|-

| 0
| <math display="inline">-1.000000</math>
| <math display="inline">-1.00000000</math>
| <math display="inline">-1.000000</math>
|-

| <math display="inline">0.5</math>
| <math display="inline">-1.224747</math>
| <math display="inline">-1.22474487</math>
| <math display="inline">-1.224747</math>
|-

| <math display="inline">1.0</math>
| <math display="inline">-1.414217</math>
| <math display="inline">-1.41421356</math>
| <math display="inline">-1.414217</math>
|-

| <math display="inline">1.5</math>
| <math display="inline">-1.581147</math>
| <math display="inline">-1.58113883</math>
| <math display="inline">-1.581147</math>
|-

| <math display="inline">2.0</math>
| <math display="inline">-1.732057</math>
| <math display="inline">-1.73205081</math>
| <math display="inline">-1.732057</math>
|}

<span id='s0025'></span>
==5. Results and discussion==

The homotopy analysis method is a powerfull technique to solve nonlinear partial and ordinary differential equations with initial and boundary conditions. The set of nonlinear partial differential equations [[#e0040|(8)]] and [[#e0045|(9)]] with boundary conditions [[#e0050|(10)]] and [[#e0055|(11)]] is solved analytically by means of homotopy analysis method and the effects of involved parameters are illustrated through graphs. In this section, we are going to present the detail analysis of flow and heat transfer of couple stress fluid over oscillating stretching surface in the presence of heat source and sink.

[[#f0020|Fig. 4]](a) and (b) are plotted to show the behavior of couple stress parameter ''K  '' and magneto-porous parameter <math display="inline">\beta </math> on the time series of velocity <math display="inline">f^{{'}}</math> at a fixed location <math display="inline">y=0.25</math> from the sheet, respectively. In [[#f0020|Fig. 4]](a) the effects of couple stress parameter ''K  '' are illustrated by keeping <math display="inline">S=1</math> and <math display="inline">\beta =0.1</math>. From this figure we observe that an increase in couple stress parameter ''K'' results in an increase in amplitude of velocity. This is perhaps due to the presence of couple stresses.  [[#f0020|Fig. 4]](b) points out the effects of magneto-porous parameter <math display="inline">\beta </math> by keeping other parameters constant. This figure shows that the amplitude of velocity decreases by increasing magneto-porous parameter ''β''.

<span id='f0020'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr4.jpg|center|355px|Time series of the velocity profile f′ at a distance y=0.25 from the sheet in ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 4.

Time series of the velocity profile <math display="inline">f^{{'}}</math> at a distance <math display="inline">y=0.25</math> from the sheet in the time period <math display="inline">\tau \in \left[0\mbox{,}10\pi \right]</math>: (a) The effects of couple stress parameter ''K  '' with <math display="inline">S=1</math> and <math display="inline">\beta =0.1</math>, (b) The effects of magneto-porous parameter <math display="inline">\beta </math> with <math display="inline">S=0.5</math> and <math display="inline">K=0.1</math>.

</span>
|}

[[#f0025|Fig. 5]](a)–(d) demonstrate the effects of couple stress parameter ''K  '' on velocity profile <math display="inline">f^{{'}}</math> at four different time instants <math display="inline">\tau =8.5\pi \mbox{,}\tau =9\pi \mbox{,}\tau =9.5\pi </math> and <math display="inline">\tau =10\pi </math>. [[#f0025|Fig. 5]](a) characterizes the effects of couple stress parameter ''K  '' on the velocity profile <math display="inline">f^{{'}}</math> at <math display="inline">\tau =8.5\pi </math>. This figure depicts that the velocity <math display="inline">f^{{'}}</math> decrease as we increase ''K''. This decrease in velocity results in decrease in boundary layer thickness. In  [[#f0025|Fig. 5]](b), the behavior of couple stress parameter ''K  '' on the <math display="inline">f^{{'}}</math> is shown at <math display="inline">\tau =9\pi </math>. This figure shows an increase in magnitude of velocity at this instant. The effects of ''K  '' at <math display="inline">\tau =9.5\pi </math> and <math display="inline">\tau =10\pi </math> are illustrated in [[#f0025|Fig. 5]](c) and (d), respectively. It is observed from these figures that at <math display="inline">\tau =9.5\pi </math> and <math display="inline">\tau =10\pi </math>, the velocity <math display="inline">f^{{'}}</math> decreases with ''K  ''. It is interesting to note that the effects of couple stress are to the increase the amplitude of velocity when <math display="inline">f^{{'}}=1</math> at the surface. For other instants when <math display="inline">f^{{'}}=-1</math> or 0, the magnitude of velocity decreases due to the presence of couple stresses.

<span id='f0025'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr5.jpg|center|px|Transverse profiles of the velocity field f′ for four different values of couple ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 5.

Transverse profiles of the velocity field <math display="inline">f^{{'}}</math> for four different values of couple stress parameter ''K  '' in the fifth period <math display="inline">\tau \in \left[8\pi \mbox{,}10\pi \right]</math> for which a periodic velocity field has been reached: (a) <math display="inline">\tau =8.5\pi </math>, (b) <math display="inline">\tau =9\pi </math>, (c) <math display="inline">\tau =9.5\pi </math> and (d) <math display="inline">\tau =10\pi </math>, with <math display="inline">S=1</math> and <math display="inline">\beta =0.1</math>.

</span>
|}

The effects of magneto-porous parameter <math display="inline">\beta </math> on velocity are illustrated in [[#f0030|Fig. 6]](a)–(d) at four different times keeping other parameters fixed. In [[#f0030|Fig. 6]](a) the effects of magneto-porous parameter <math display="inline">\beta </math> are plotted at <math display="inline">\tau =8.5\pi </math>. From this figure it can be examined that velocity decreases with increase in magneto-porous parameter <math display="inline">\beta </math>. The effects of magneto-porous parameter on velocity at other three time instants are illustrated in [[#f0030|Fig. 6]](c) and (d). Here, again it is observed that the magnitude of velocity decreases by increasing <math display="inline">\beta </math>.

<span id='f0030'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr6.jpg|center|px|Transverse profiles of the velocity field f′ for four different values of ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 6.

Transverse profiles of the velocity field <math display="inline">f^{{'}}</math> for four different values of magneto-porous parameter <math display="inline">\beta </math> in the fifth period <math display="inline">\tau \in \left[8\pi \mbox{,}10\pi \right]</math> for which a periodic velocity field has been reached: (a) <math display="inline">\tau =8.5\pi </math>, (b) <math display="inline">\tau =9\pi </math>, (c) <math display="inline">\tau =9.5\pi </math>, and (d) <math display="inline">\tau =10\pi </math>, with <math display="inline">S=0.5</math> and <math display="inline">K=0.1</math>.

</span>
|}

The variation of skin friction coefficient with couple stress parameter ''K  '' and magneto-porous parameter <math display="inline">\beta </math> in the first five periods <math display="inline">\tau \in [0\mbox{,}10\pi ]</math> is shown in [[#f0035|Fig. 7]](a) and (b), respectively. [[#f0035|Fig. 7]](a) shows the time series of skin friction coefficient for different values of couple stress parameter ''K  '' by taking <math display="inline">S=0.1</math> and <math display="inline">\beta =0.1</math>. From this figure we observe that skin friction coefficient varies periodically due to periodic motion of the surface and its amplitude increases with couple stress parameter ''K  ''. This observation reflects that skin friction for a Newtonian fluid flowing over an oscillatory stretching sheet is less in comparison with its value for a couple stress fluid performing the same motion. The effects of magneto-porous parameter <math display="inline">\beta </math> on skin friction coefficient are illustrated in [[#f0035|Fig. 7]](b). This figure depicts that amplitude of skin friction coefficient increases by increasing magneto-porous parameter <math display="inline">\beta </math>. Such effects of <math display="inline">\beta </math> and ''M'' are expected because the fact that magnetic force and drag force offered by the porous medium acts as a resistance to the flow.

<span id='f0035'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr7.jpg|center|355px|Time series of the skin-friction coefficient Rex1/2Cf in the first five periods ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 7.

Time series of the skin-friction coefficient <math display="inline">{Re}_x^{1/2}C_f</math> in the first five periods <math display="inline">\tau \in [0\mbox{,}10\pi ]</math>: (a) effects of couple stress parameter ''K  '' with <math display="inline">S=0.1</math> and <math display="inline">\beta =0.1</math> and (b) effects of magneto-porous parameter <math display="inline">\beta </math> with <math display="inline">S=0.5</math> and <math display="inline">K=0.1</math>.

</span>
|}

[[#f0040|Figure 8]], [[#f0045|Figure 9]], [[#f0050|Figure 10]], [[#f0055|Figure 11]] and [[#f0060|Figure 12]] present the effects of Prandtl number ''Pr  '', magneto-porous parameter <math display="inline">\beta </math>, couple stress parameter ''K  '' and heat source/sink parameter <math display="inline">\lambda </math> on the temperature field. [[#f0040|Fig. 8]] reveals the effects of Prandtl number ''Pr'' on temperature profile for both heat source and sink cases. It is seen that temperature is a decreasing function of ''Pr'' in both cases. The fact is that the increase in Prandtl number results in low thermal conductivity, as a result conduction as well as thermal boundary layer thickness decreases and hence we observe a decrease in temperature. It is observed that the presence of heat sink enhances the decrease in temperature by increasing Prandtl number. In  [[#f0045|Fig. 9]], the influence of magneto-porous parameter <math display="inline">\beta </math> on temperature field is explained graphically by keeping other parameters fixed. This figure shows that temperature increases with magneto-porous parameter <math display="inline">\beta </math>. The influence of couple stress parameter ''K'' on temperature is shown in  [[#f0050|Fig. 10]], which depicts that increase in couple stress parameter results in increase in temperature. Since we are discussing the temperature field in the presence of heat source/sink so effects of this parameter are very important for both the cases. The effects of heat source/sink parameter are illustrated in  [[#f0055|Figure 11]] and [[#f0060|Figure 12]] by keeping <math display="inline">Pr=0.4\mbox{,}\beta =0.5\mbox{,}K=0.1\mbox{,}S=0.1</math>. It is noted from [[#f0055|Fig. 11]] that as we increase the strength of the heat source, the temperature increases. This is due to the fact that heat source can add more heat to the stretching sheet which increases its temperature and thus the temperature of fluid rises. Moreover, the thermal boundary layer is also found to increase by increasing the strength of heat source. [[#f0060|Fig. 12]] is plotted to observe the effects of heat sink by keeping other parameters constant. This figure shows opposite results i.e., the temperature decreases by increasing the strength of heat sink. This is because, when strength of heat sink increases, more heat is removed from the sheet and thus a decrease in thermal boundary-layer thickness and temperature is observed. This result is of key importance for the flows where heat transfer is of prime importance. The numerical values of Nusselt number are shown in [[#t0015|Table 3]]. It is found that the Nusselt number increases with Prandtl number while it decreases with increase in heat source, couple stress and magneto-porous parameter.

<span id='f0040'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr8.jpg|center|px|The effects of Prandtl number Pr on temperature field θ at τ=8.5π with ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 8.

The effects of Prandtl number ''Pr  '' on temperature field <math display="inline">\theta </math> at <math display="inline">\tau =8.5\pi </math> with <math display="inline">K=0.1\mbox{,}S=0.1\mbox{,}\beta =0.2</math>, (a) in heat source case <math display="inline">\lambda =0.3</math>, (b) in heat sink case <math display="inline">\lambda =-0.3</math>.

</span>
|}

<span id='f0045'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr9.jpg|center|311px|The effects of magneto-porous parameter β on temperature field θ at τ=8.5π with ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 9.

The effects of magneto-porous parameter <math display="inline">\beta </math> on temperature field <math display="inline">\theta </math> at <math display="inline">\tau =8.5\pi </math> with <math display="inline">Pr=1.5\mbox{,}K=0.1\mbox{,}S=0.1\mbox{,}\lambda =0.3</math>.

</span>
|}

<span id='f0050'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr10.jpg|center|311px|The effects of couple stress number K on temperature field θ at τ=8.5π with ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 10.

The effects of couple stress number ''K  '' on temperature field <math display="inline">\theta </math> at <math display="inline">\tau =8.5\pi </math> with <math display="inline">Pr=1.5\mbox{,}\beta =5\mbox{,}\lambda =0.4\mbox{,}S=</math><math>0.6</math>.

</span>
|}

<span id='f0055'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr11.jpg|center|311px|The effects of heat source parameter λ on temperature field θ at τ=8.5π with ...]]


|-
| <span style="text-align: center; font-size: 75%;">

Figure 11.

The effects of heat source parameter <math display="inline">\lambda </math> on temperature field <math display="inline">\theta </math> at <math display="inline">\tau =8.5\pi </math> with <math display="inline">Pr=0.4\mbox{,}\beta =0.5\mbox{,}K=0.1\mbox{,}S=0.1</math>.

</span>
|}

<span id='f0060'></span>

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; max-width: 100%;" 
|-
|


[[Image:draft_García_824576488-1-s2.0-S1110016816000715-gr12.jpg|center|309px|The effects of heat sink parameter λ on temperature field θ at τ=8.5π with ...]]

Figure 12.

The effects of heat sink parameter <math display="inline">\lambda </math> on temperature field <math display="inline">\theta </math> at <math display="inline">\tau =8.5\pi </math> with <math display="inline">Pr=0.4\mbox{,}\beta =0.5\mbox{,}K=0.1\mbox{,}S=0.1</math>.


|-
| 
|}

<span id='t0015'></span>

{| class="wikitable" style="min-width: 60%;margin-left: auto; margin-right: auto;"
|+

Table 3.

The numerical values of the Local Nusselt number for different values of Prandtl number ''Pr'', couple stress parameter ''K  '', magneto-porous parameter <math display="inline">\beta </math> and heat source parameter <math display="inline">\lambda </math> at time <math display="inline">\tau =\pi /2</math> with <math display="inline">S=1</math>.

|-

! ''Pr''
! <math display="inline">\lambda </math>
! ''K''
! <math display="inline">\beta </math>
! <math display="inline">{Re}_x^{-1/2}{Nu}_x</math>
|-

| <math display="inline">0.5</math>
| <math display="inline">0.2</math>
| <math display="inline">5.0</math>
| <math display="inline">5.0</math>
| <math display="inline">0.596791</math>
|-

| <math display="inline">1.0</math>
| 
| 
| 
| <math display="inline">0.634885</math>
|-

| <math display="inline">1.5</math>
| 
| 
| 
| <math display="inline">0.671806</math>
|-

| <math display="inline">2.0</math>
| 
| 
| 
| <math display="inline">0.705075</math>
|-

| <math display="inline">0.5</math>
| <math display="inline">0.2</math>
| 
| 
| <math display="inline">0.596791</math>
|-

| 
| <math display="inline">0.4</math>
| 
| 
| <math display="inline">0.547204</math>
|-

| 
| <math display="inline">0.6</math>
| 
| 
| <math display="inline">0.496478</math>
|-

| 
| <math display="inline">0.8</math>
| 
| 
| <math display="inline">0.44459</math>
|-

| 
| <math display="inline">0.2</math>
| <math display="inline">0.2</math>
| 
| <math display="inline">0.573683</math>
|-

| 
| 
| <math display="inline">0.4</math>
| 
| <math display="inline">0.573021</math>
|-

| 
| 
| <math display="inline">0.8</math>
| 
| <math display="inline">0.572122</math>
|-

| 
| 
| <math display="inline">1.0</math>
| 
| <math display="inline">0.571885</math>
|-

| 
| 
| <math display="inline">5.0</math>
| <math display="inline">0.5</math>
| <math display="inline">0.643221</math>
|-

| 
| 
| 
| <math display="inline">1.5</math>
| <math display="inline">0.630873</math>
|-

| 
| 
| 
| <math display="inline">2.5</math>
| <math display="inline">0.619684</math>
|-

| 
| 
| 
| <math display="inline">3.5</math>
| <math display="inline">0.609657</math>
|}

<span id='s0030'></span>
==6. Concluding remarks==

Flow and heat transfer of a couple stress fluid over an oscillatory stretching sheet in the presence of heat source/sink is investigated using HAM. It is observed that flow and heat transfer characteristics are greatly influenced by the presence of couple stress and heat source/sink parameter. In fact, the skin friction coefficient increases by increasing couple stress and magneto-porous parameter. It is also interesting to note that rate of heat transfer is enhanced by increasing magneto-porous parameter, heat source and couple stress parameter, while it decreases for large values of Prandtl number and heat sink parameter. Such observations may have interesting implications where it is desired to reduce skin friction and enhance heat transfer.

==Acknowledgments==

We are grateful to the reviewers for their useful comments. The second author is grateful to the Higher Education Commission of Pakistan for financial assistance.

==References==

<ol style='list-style-type: none;margin-left: 0px;'><li><span id='b0005'></span>
[[#b0005|[1]]] B.C. Sakiadis; Boundary layer behavior on continuous solid surfaces: I boundary layer equations for two dimensional and axisymmetric flow; AIChE J., 7 (1961), pp. 26–28</li>
<li><span id='b0010'></span>
[[#b0010|[2]]] L.J. Crane; Flow past a stretching plate; Z. Angew. Math. Phy. (ZAMP), 21 (1970), pp. 645–647</li>
<li><span id='b0015'></span>
[[#b0015|[3]]] K.R. Rajagopal, T.Y. Na, A.S. Gupta; Flow of a visco-elastic fluid over a stretching sheet; Rheol. Acta, 23 (1984), pp. 213–215</li>
<li><span id='b0020'></span>
[[#b0020|[4]]] P.S. Gupta, A.S. Gupta; Heat and mass transfer on a stretching sheet with suction or blowing; Can. J. Chem. Eng., 55 (1977), pp. 744–746</li>
<li><span id='b0025'></span>
[[#b0025|[5]]] B.S. Dandapat, A.S. Gupta; Flow and heat transfer in a viscoelastic fluid over a stretching sheet; Int. J. Nonlinear Mech., 24 (1989), pp. 215–219</li>
<li><span id='b0030'></span>
[[#b0030|[6]]] J.B. McLeod, K.R. Rajagopal; On the uniqueness of flow of a Navier-Stokes fluid due to stretching boundries; Arch. Rat. Mech. Anal., 98 (1987), p. 385</li>
<li><span id='b0035'></span>
[[#b0035|[7]]] D. Rollins, K. Vajravelu; Heat transfer in a second-order fluid over a continuous stretching surface; Acta Mech., 89 (1991), pp. 167–178</li>
<li><span id='b0040'></span>
[[#b0040|[8]]] K.R. Rajagopal, R.K. Bhatnagar; Exact solutions for some simple flows of an Oldroyd-B fluid; Acta Mech., 113 (1995), pp. 233–239</li>
<li><span id='b0045'></span>
[[#b0045|[9]]] R. Cortell; A note on flow and heat transfer of a viscoelastic liquid over a stretching sheet; Int. J. Non-Linear Mech., 41 (2006), pp. 78–85</li>
<li><span id='b0050'></span>
[[#b0050|[10]]] R. Cortell; Flow and heat transfer of an electrically conducting fluid of second grade over a stretching sheet subject to suction and to a transverse magnetic field; Int. J. Heat Mass Transfer, 49 (2006), pp. 1851–1856</li>
<li><span id='b0055'></span>
[[#b0055|[11]]] R. Nazar, N. Amin, D. Filip, I. Pop; Stagnation point flow of a micropolar fluid towards a stretching sheet; Int. J. Non-Linear Mech., 39 (2004), pp. 1227–1235</li>
<li><span id='b0060'></span>
[[#b0060|[12]]] A. Ishak, R. Nazar, I. Pop; Mixed convection stagnation point flow of a micropolar fluid towards a stretching sheet; Mecccanica, 43 (2008), pp. 411–418</li>
<li><span id='b0065'></span>
[[#b0065|[13]]] P.D. Ariel; Axisymmetric flow of a second grade fluid past a stretching sheet; Int. J. Eng. Sci., 39 (2001), pp. 529–553</li>
<li><span id='b0070'></span>
[[#b0070|[14]]] T. Hayat, M. Sajid; Analytic solution for axisymmetric flow and heat transfer of a second grade fluid past a stretching sheet; Int. J. Heat Mass Transfer, 50 (2007), pp. 75–84</li>
<li><span id='b0075'></span>
[[#b0075|[15]]] T. Hayat, Z. Abbas, M. Sajid; Heat and mass transfer analysis on the flow of a second grade fluid in the presence of chemical reaction; Phys. Lett. A, 372 (2008), pp. 2400–2408</li>
<li><span id='b0080'></span>
[[#b0080|[16]]] T. Hayat, T. Javed, Z. Abbas; Slip flow and heat transfer of a second grade fluid past a stretching sheet through a porous space; Int J. Heat Mass Transfer, 51 (2008), pp. 4528–4534</li>
<li><span id='b0085'></span>
[[#b0085|[17]]] M. Sajid, I. Ahmad, T. Hayat, M. Ayub; Unsteady flow and heat transfer of a second grade fluid over a stretching sheet; Commun. Nonlinear Sci. Numer. Simul., 14 (2009), pp. 96–108</li>
<li><span id='b0090'></span>
[[#b0090|[18]]] T. Hayat, T. Javed, M. Sajid; Analytic solution for rotating flow and heat transfer analysis of a third-grade fluid; Acta Mech., 191 (2007), pp. 219–229</li>
<li><span id='b0095'></span>
[[#b0095|[19]]] M. Sajid, T. Hayat; Non-similar series solution for boundary layer flow of a third-order fluid over a stretching sheet; Appl. Math. Comput., 189 (2007), pp. 1576–1585</li>
<li><span id='b0100'></span>
[[#b0100|[20]]] T. Hayat, M. Mustafa, I. Pop; Heat and mass transfer for Soret and Dufour’s effect on mixed convection boundary layer flow over a stretching vertical surface in a porous medium filled with a viscoelastic fluid; Commun. Nonlinear Sci. Numer. Simul., 15 (2010), pp. 1183–1196</li>
<li><span id='b0105'></span>
[[#b0105|[21]]] M. Turkyilmazoglu; The analytical solution of mixed convection heat transfer and fluid flow of a MHD viscoelastic fluid over a permeable stretching surface; Int. J. Mech. Sci., 77 (2013), pp. 263–268</li>
<li><span id='b0110'></span>
[[#b0110|[22]]] M. Turkyilmazoglu; Three dimensional MHD flow and heat transfer over a stretching/shrinking surface in a viscoelastic fluid with various physical effects; Int. J. Heat Mass Transfer, 78 (2014), pp. 150–155</li>
<li><span id='b0115'></span>
[[#b0115|[23]]] M. Turkyilmazoglu; Multiple solutions of hydromagnetic permeable flow and heat for viscoelastic fluid; J. Thermophys. Heat Transfer, 25 (2011), pp. 595–605</li>
<li><span id='b0120'></span>
[[#b0120|[24]]] M. Turkyilmazoglu; Analytic heat and mass transfer of the mixed hydrodynamic/thermal slip MHD viscous flow over a stretching sheet; Int. J. Mech. Sci., 53 (2011), pp. 886–896</li>
<li><span id='b0125'></span>
[[#b0125|[25]]] M. Turkyilmazoglu; Multiple solutions of heat and mass transfer of MHD slip flow for the viscoelastic fluid over a stretching sheet; Int. J. Therm. Sci., 50 (2011), pp. 2264–2276</li>
<li><span id='b0130'></span>
[[#b0130|[26]]] M. Turkyilmazoglu; An analytical treatment for the exact solutions of MHD flow and heat over two–three dimensional deforming bodies; Int. J. Heat Mass Transfer, 90 (2015), pp. 781–789</li>
<li><span id='b0135'></span>
[[#b0135|[27]]] C.Y. Wang; Nonlinear streaming due to the oscillatory stretching of a sheet in a viscous; Acta Mech., 72 (1988), pp. 261–268</li>
<li><span id='b0140'></span>
[[#b0140|[28]]] B. Siddappa, S. Abel, V. Hongunti; Oscillatory motion of a viscoelastic fluid past a stretching sheet; ll Nuovo Cimento D, 17 (1995), pp. 53–60</li>
<li><span id='b0145'></span>
[[#b0145|[29]]] Z. Abbas, Y. Wang, T. Hayat, M. Oberlack; Slip effects and heat transfer analysis in a viscous fluid over an oscillatory stretching surface; Int. J. Numer. Methods Fluids, 59 (2009), pp. 443–458</li>
<li><span id='b0150'></span>
[[#b0150|[30]]] Z. Abbas, Y. Wang, T. Hayat, M. Oberlack; Hydromagnetic flow in a viscoelastic fluid due to the oscillatory stretching surface; Int. J. Nonlinear Mech., 43 (2008), pp. 783–797</li>
<li><span id='b0155'></span>
[[#b0155|[31]]] L.C. Zheng, X. Jin, X.X. Zhang, J.H. Zhang; Unsteady heat and mass transfer in MHD flow over an oscillatory stretching surface with Soret and Dufour effects; Acta Mech. Sinica, 29 (5) (2013), pp. 667–675</li>
<li><span id='b0160'></span>
[[#b0160|[32]]] N. Ali, S.U. Khan, Z. Abbas; Hydromagnetic flow and heat transfer of a jeffrey fluid over an oscillatory stretching surface; Z. Naturforsch. A, 70 (7) (2015), pp. 567–576</li>
<li><span id='b0165'></span>
[[#b0165|[33]]] S.U. Khan, N. Ali, Z. Abbas; Hydromagnetic flow and heat transfer over a porous oscillating stretching surface in a viscoelastic fluid with porous medium; PLOS ONE (2015) [http://dx.doi.org/10.1371/journal.pone.0144299 http://dx.doi.org/10.1371/journal.pone.0144299] December 14</li>
<li><span id='b0170'></span>
[[#b0170|[34]]] V.K. Stokes; Couple Stress Fluids, 9 (1966), pp. 1709–1715</li>
<li><span id='b0175'></span>
[[#b0175|[35]]] D. Srinivasacharya, K. Kaladhar; Soret and Dufour effects in a mixed convection couple stress fluid with heat and mass fluxes; Lat. Am. Appl. Res., 41 (2011), pp. 353–358</li>
<li><span id='b0180'></span>
[[#b0180|[36]]] L.M. Srivastava; Peristaltic transport of a couple stress fluid; Rheol. Acta, 25 (1986), pp. 638–641</li>
<li><span id='b0185'></span>
[[#b0185|[37]]] E.F. El-Shehawey, K.S. Mekheimer; Couple-stresses in peristaltic transport of fluids; J. Phys. D: Appl. Phys., 27 (1994), p. 1163</li>
<li><span id='b0190'></span>
[[#b0190|[38]]] D. Pal, N. Rudraiah, R. Devanathan; A couple stress model of blood flow in the microcirculation; Bull. Math. Biol., 50 (1988), p. 329</li>
<li><span id='b0195'></span>
[[#b0195|[39]]] H.L. Chiang, J.R. Lin, C.H. Hsu, Y.P. Chang; Linear stability analysis of a rough short journal bearing lubricated with non-newtonian fluids; Tribol. Int., 17 (2004), pp. 867–877</li>
<li><span id='b0200'></span>
[[#b0200|[40]]] N.B. Naduvinamani, P.S. Hiremath, G. Gurubasavaraj; Effect of surface roughness on the couple-stress squeeze film between a sphere and a flat plate; Tribol. Int., 38 (2005), pp. 451–458</li>
<li><span id='b0205'></span>
[[#b0205|[41]]] C.W.C. Jian, C.K. Chen; Chaos and bifurcation of a flexible rotor supported by porous squeeze couple stress fluid film journal bearings with non-linear suspension; Chaos Solit. Fract., 35 (2008), pp. 358–375</li>
<li><span id='b0210'></span>
[[#b0210|[42]]] R.F. Lu, J.R. Lin; A theoretical study of combined effects of non-Newtonian rheology and viscosity-pressure dependence in the sphere-plate squeeze-film system; Tribo. In., 40 (1) (2007), pp. 125–131</li>
<li><span id='b0215'></span>
[[#b0215|[43]]] S.J. Liao; Advance in the Homotopy Analysis Method, vol. 1, World Scientific Publishing, 5 Toh Tuck Link, Singapore (2014)</li>
<li><span id='b0220'></span>
[[#b0220|[44]]] S. Abbasbandy; The application of homotopy analysis method to nonlinear equations arising in heat transfer; Phys. Lett. A, 360 (2006), pp. 109–113</li>
<li><span id='b0225'></span>
[[#b0225|[45]]] M. Turkyilmazoglu; Analytic approximate solutions of rotating disk boundary layer flow subject to a uniform suction or injection; Int. J. Mech. Sci., 52 (2010), pp. 1735–1744</li>
<li><span id='b0230'></span>
[[#b0230|[46]]] S. Abbasbandy; Homotopy analysis method for heat radiation equations; Int. Commun. Heat Mass Transfer, 34 (2007), pp. 380–387</li>
<li><span id='b0235'></span>
[[#b0235|[47]]] I.A. Abdallah; Homotopy analytical solution of MHD fluid flow and heat transfer problem; Appl. Math. Inform. Sci., 3 (2009), pp. 223–233</li>
<li><span id='b0240'></span>
[[#b0240|[48]]] T. Hayat, M. Sajid; Homotopy analysis of MHD boundary layer flow of an upper-convected Maxwell fluid; Int. J. Eng. Sci., 45 (2007), pp. 393–401</li>
<li><span id='b0245'></span>
[[#b0245|[49]]] B. Raftari, K. Vajravelu; Homotopy analysis method for MHD viscoelastic fluid flow and heat transfer in a channel with a stretching wall; Commun. Nonlinear Sci. Numer. Simul., 17 (2012), p. 4149</li>
<li><span id='b0250'></span>
[[#b0250|[50]]] S.J. Liao; A uniformly valid analytic solution of 2D viscous flow past a semi-infinite flat plate; J. Fluid Mech., 385 (1999), pp. 101–128</li>
<li><span id='b0255'></span>
[[#b0255|[51]]] S.J. Liao; An explicit, totally analytic approximation of Blasius’ viscous flow problems totally analytic approximation of Blasius’ viscous flow problems; Int. J. Non-Linear Mech., 34 (4) (1999), pp. 759–778</li>
<li><span id='b0260'></span>
[[#b0260|[52]]] S.J. Liao; On the analytic solution of magnetohydrodynamic flows of non-Newtonian fluids over a stretching sheet; J. Fluid Mech., 488 (2003), pp. 189–212</li>
<li><span id='b0265'></span>
[[#b0265|[53]]] M. Turkyilmazoglu; Solution of the Thomas-Fermi equation with a convergent approach; Commun. Nonlinear Sci. Numer. Simul., 17 (2012), pp. 4097–4103</li>
</ol>