\documentclass[reqno]{amsart}
\usepackage{hyperref}
\usepackage{graphicx}

\AtBeginDocument{{\noindent\small
International Conference on Applications of Mathematics to Nonlinear Sciences,\newline
\emph{Electronic Journal of Differential Equations},
Conference 24 (2017), pp. 85--101.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu}
\thanks{\copyright 2017 Texas State University.}
\vspace{8mm}}

\begin{document} \setcounter{page}{85}
\title[\hfilneg EJDE-2017/conf/24\hfil 
 Free surface dynamics]
{Free surface dynamics of thin MHD second-grade fluid over a
heated nonlinear stretching sheet}

\author[K. K. Patra, S. Panda, M. Sellier \hfil EJDE-2017/conf/24\hfilneg]
{Kiran Kumar Patra, Satyananda Panda, Mathieu Sellier}

\address{Kiran Kumar Patra \newline
Department of Mathematics,
National Institute of Technology Calicut,
NIT(P.O)-673601, Kerala, India}
\email{kirankumarpatra1984@gmail.com}

\address{Satyananda Panda (corresponding author)\newline
Department of Mathematics,
National Institute of Technology Calicut,
NIT(P.O)-673601, Kerala, India}
\email{satyanand@nitc.ac.in }

\address{Mathieu Sellier \newline
Department of Mechanical Engineering,
University of Canterbury, Private Bag 4800,
 Christchurch 8140, New Zealand}
\email{mathieu.sellier@canterbury.ac.nz}

\thanks{Published November 15, 2017.}
\subjclass[2010]{76A05, 76A10, 76A20, 76M12, 80A20}
\keywords{Thin liquid film; heat transfer; second-grade fluid;
 free surface flow; 
\hfill\break\indent magnetic field; long-wave theory}

\begin{abstract}
 This article presents a long-wave theory for the free surface dynamics
 of magnetohydrodynamics (MHD) second-grade fluid over a non-uniform heated
 flat elastic sheet. An evolution equation for the film thickness is derived
 from the instationary Navier-Stokes equations using regular asymptotic
 expansion with respect to the small aspect ratio of the flow domain.
 The derived thin film equation is solved numerically using finite volume
 method on a uniform grid system with implicit flux discretization.
 The finding reveals the dependency of the thinning behavior of the fluid
 film on the stretching speed and the non-Newtonian second-grade parameter.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\allowdisplaybreaks

\section{Introduction}\label{sec:sect1}

 Many theoretical studies on the flow of fluid and heat transfer due to stretching
sheet in the context of polymer extrusion, continuous casting, drawing of
the plastic sheets, cable coating, etc., can be found in the literature, about, e.g.
Newtonian fluid
 \cite{Andersson2000, Aziz2010, Liu2008,dandapat,vajravelu,  Khan2011,salleh, Wang},
non-Newtonian fluid 
\cite{andersson1996flow, hayat2007analytic,sajid2009unsteady, sajid2008influence} and
 MHD effect \cite{hayat2008influence,ali,Noor2010} as well as derivation
of boundary layer equations. In general, such model reductions are based
on the uniform film thickness assumption
which enables the similarity transformation. Thereby the set of partial
differential equations are reduced to a more tractable one of ordinary
differential equations.

Recognizing the restrictions of the plane interface assumptions,
Dandapat and co-workers were the first to
extend the formulation to account for local deformation of the free surface
in \cite{dandapatZAMP2006, dandapatPoF2006}. The authors exploit the slenderness
of the flow domain to derive a long-wave approximation of
the Navier-Stokes equations and solve the resulting governing
equation using the matched asymptotic method. Lately, this
work was extended to include the heat transfer problem
\cite{ dandapat2011thin,santra2009}.
In this work, the nonuniform temperature distribution at
the stretching sheet induces an inhomogeneous temperature
field in the film. Consequently, a surface temperature gradient
develops at the film free surface. As a result of the surface
tension gradients, the film thickness varies along the flow, and
these deformations are advected in the stretching direction. Also,
a free surface model based on a long-wave theory for the thin film
dynamics of Casson fluid over a nonlinear
stretching sheet including magnetic effect has been
recently deduced in \cite{singh}.


This work focuses on the systematic derivation, in the spirit of
\cite{dandapat2011thin}, of the thin film equation for
a second-grade non-Newtonian MHD fluid over a heated steady stretching sheet without
the restriction of the plane interface assumption. One motivation for this study
is the flow of mucus in biological tissues which undergo expansion or contraction.
A particular example is pulmonary alveoli which are
covered with a lining of non-Newtonian fluid \cite{levy} and which undergo periodic
expansion and contraction.

This article is organized as follows. The mathematical model for the flow of
second-grade fluid is described in Sec.~\eqref{sec:sect2}.
The long-wave theory by using the standard expansion technique with
respect to a small aspect ratio of the
flow domain for model reduction is presented in Sec.~\eqref{sec:sect3}. The
thin film model equation is given in Sec.~\eqref{sec:sect4}.
 The numerical procedure for the numerical solution of the thin
film equation is explained in Sec.~\eqref{sec:sect5}. 
In Sec.~\eqref{sec:sect6} we
discuss the numerical investigation of the thin film
equation, and give some conclusions from this work.

 \section{Physical model and problem formulation}\label{sec:sect2}

We consider here an unsteady planar fluid which lies over a heated
stretching sheet in the presence
of transverse magnetic field $\mathbf{B}_0$ as shown in
Figure \eqref{fig:flow_geometry}.
 The elastic sheet lies at $z=0$, and
the liquid-gas interface lies at $z=h(x,t)$, where the $x$-axis is directed
along the stretching layer, and the
$z$-axis is normal to the sheet in the outward direction toward
the fluid. Gravity acts along the negative z-direction.
 Further, the
surface at $z$ = 0 starts stretching from rest and within a very
short time attains the stretching velocity $u = U(x)$.
The evaporation and buoyancy effects are neglected considering
the liquid is non-volatile and thin.
The elastic surface is heated with non-uniform temperature $T_{s}$ a function of
$x$ alone, and the ambient gas phase is at constant temperature $T_{a} $.
A constant uniform magnetic field of strength $B_0$ is applied transversely
in the parallel direction to the $z$-axis. Since the sheet is nonconducting,
an uniform electric field is applied along
the z-direction to generate Lorentz force.

The fluid is assumed to be incompressible and non-Newtonian second-grade fluid with
constant viscosity $\mu$, density $\rho$, specific heat $c_{p}$, thermal
conductivity $k$ and the electric conductivity $\lambda$. The surface tension
of the liquid-gas interface decreases linearly with temperature $T$ according to
\begin{equation}\label{eq:st}
\sigma = \sigma_{a}\big(1-\gamma(T - T_{a})\big),
\end{equation}
where $\sigma_{a}$ is the surface tension at $T=T_{a}$, and
$\gamma=-\frac{1}{\sigma_a}\big(\frac{d\sigma}{dT}\big)_{T=T_a}$
a positive constant specific to the fluid.

\begin{figure}[hbt]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig1} % pd.eps Prob_desc}
\end{center}
\caption{Sketch of flow geometry}
\label{fig:flow_geometry}
\end{figure}

Neglecting viscous dissipation and radiation effects,
the motion of the second-grade fluid due to the heated stretching sheet is
governed by the continuity equation of mass flow, the momentum equation
and the temperature equation.

 \begin{itemize}
 \item Continuity equation:
 \begin{equation}\label{eq:contd}
 \nabla \cdot \mathbf{V} = \mathbf{0};
 \end{equation}

 \item Momentum equation:
 \begin{equation}\label{eq:momd}
 \rho \frac{D\mathbf{V}}{Dt} = \nabla \cdot {\boldsymbol{\tau}}
+ \rho \mathbf{g} + \mathbf{F};
 \end{equation}

 \item Temperature equation:
 \begin{equation}\label{eq:tempd}
 \rho c_p \frac{DT}{Dt}=k\nabla^2T;
 \end{equation}

\item Free surface boundary conditions i.e. at $z=h(x,t)$:
\begin{gather}\label{eq:bc2d}
p_{a} + \hat{\mathbf{n}}\cdot \boldsymbol{\tau}\cdot \hat{\mathbf{n}}
 = -\sigma(\boldsymbol{\nabla}\cdot \hat{\mathbf{n}}), \\
\label{eq:bc3d}
\hat{\mathbf{n}}\cdot \boldsymbol{\tau} \cdot \hat{\mathbf{t}}
= \nabla \sigma \cdot \hat{\mathbf{t}}, \\
\label{eq:bc4d}
h_{t} + u h_{x} = w, \\
q_i = -k \nabla T \cdot \hat{\textbf{n}} = \alpha\left(T-T_{a}\right)
\end{gather}

\item Boundary conditions on the stretching sheet i.e. at $z=0$:
\begin{equation}\label{eq:bc1d}
u(x, 0, t) = U(x), \quad  w(x, 0, t) = 0, \quad T(x, 0, t) = T_{s}(x);
\end{equation}

\item Initial condition:
\begin{equation}\label{eq:ic}
u(x,z,0) = w(x, z,0) = 0, \quad h(x,0) = h_0 +\delta(x)
\end{equation}
 \end{itemize}

The symbols $\mathbf{V}(x,z,t)=(u(x,z,t), w(x,z,t))$,
${\boldsymbol{\tau}}$, $\mathbf{g}$ and $\mathbf{F}$
 denote the fluid velocity at position $(x,z)$ and
time $t$, the Cauchy stress tensor, the acceleration due to gravity and the
Lorentz force.
The free surface boundary conditions are due to balance of stresses,
the kinematic condition and the convective heat flux at the interface.
Here, $\hat{\mathbf{n}}$ and $\hat{\mathbf{t}}$ are the unit normal and
tangential vectors on the surface, respectively. The $p_{a}$ stands for the
atmospheric pressure at the free surface, $\sigma$ is the surface tension of the
fluid, $q_i$ is the heat flux, $\alpha$ is the rate of heat transfer from
the liquid to the ambient gas phase and the symbol $\nabla$ stands for the
gradient operator. The subscripts $x$ and $t$ stand for the partial
differentiation with respect to $x$ and $t$ respectively. The symbol
$h_0$ is the characteristic height of the free surface and $\delta$
is the initial small disturbance from $h_0$.

The \emph{Cauchy-stress tensor}, given by Rivlin and Ericksen \cite{rivlin}
for a second-grade fluid can be written as
\begin{equation}\label{eq:cst}
\mathbf{\boldsymbol{\tau}} = - p\mathbf{I} + \mu \mathbf{A}_1 
+ \alpha_1 \mathbf{A}_2
+ \alpha_{2} \mathbf{A}_1^{2},
\end{equation}
where $p$ is the pressure, $\mathbf{I}$ is the identity tensor.
The material constants $\alpha_1$ and $\alpha_{2}$ are the first and
second normal stress coefficients.
 The quantities $\mathbf{A}_1$ and $\mathbf{A}_2$ are the first two
 Rivlin-Ericksen tensors and they are defined as
\begin{equation} %\label{eq:r-e-tensor}
\mathbf{A}_1 = \left(\nabla \mathbf{V}\right)
 + \left(\nabla \mathbf{V}\right)^{*}, \quad
\mathbf{A}_2 = \frac{D}{Dt} \mathbf{A}_1
 + \mathbf{A}_1 \cdot \left(\nabla \mathbf{V}\right) +
\left({\nabla \mathbf{V}}\right)^{*} \cdot \mathbf{A}_1\,,
\label{eq:an}
\end{equation}
where $D/Dt$ is the material time derivative, and the superscript $(*)$
is used for the transpose.

The constitutive model \eqref{eq:cst} is derived by considering second order
approximation of retardation parameter. Dunn and Fosdick~\cite{dunn}
have shown that this model equation is invariant under transformation
and therefore the material constants must meet the following
restriction
\begin{equation}\label{eq:cstres}
\mu \ge 0, \quad \alpha_1 \ge 0, \quad \alpha_1 + \alpha_{2} = 0.
\end{equation}
Fluids characterized by these restrictions \eqref{eq:cstres} are
called second-grade fluid. The fluid model represented by \eqref{eq:cst}
with the relationship  \eqref{eq:cstres} is compatible with the hydrodynamics.
The third relations of \eqref{eq:cstres} is the
consequence of satisfying the Clausis-Duhem inequality by fluid motion and
a second relation arises due to
the assumptions that specific Helmholtz free energy of the fluid takes its
minimum value in equilibrium.
The fluid satisfying model \eqref{eq:cst} with $\alpha_i < 0 $;
 ($\alpha_i = 1, 2$) is termed as
second-order fluid and with $\alpha_i > 0$ is termed as second-grade fluid.
Although second-order fluid is obeying model
\eqref{eq:cst} with $\alpha_1< \alpha_{2}, \alpha_1 <0$,
exhibits some undesirable instability characteristic
(Fosdick and Rajagopal~\cite{fosdicketal}).
 The second order approximation is valid at low shear rate
(Dunn and Rajagopal~\cite{dunnetal}).

The fluid is flowing under the environment of uniform transverse magnetic field,
therefore momentum of the  fluid is influenced by the Lorentz force
$\mathbf{F}=\lambda(\textbf{E}+\mathbf{V}\times\mathbf{B}_0)\times\mathbf{B}_0$,
 where $\textbf{E}+\mathbf{V}\times\mathbf{B}_0$ represents the total current
density with magnetic Reynolds number
 $R_m<<1$. The electric field $\textbf{E}=\mathbf{0}$, as there is no
electric current present during the fluid flow.
 The term $\mathbf{V}\times\mathbf{B}_0$ is the potential difference across
the fluid.
 Using basic vector calculus the expression for Lorentz force can be
calculated as follows
 \begin{equation}\label{eq:efd}
 \mathbf{F}=\lambda(\mathbf{V}\times\mathbf{B}_0)\times\mathbf{B}_0
= -\lambda[\mathbf{V}(\mathbf{B}_0\cdot\mathbf{B}_0)
 -\mathbf{B}_0(\mathbf{V}\cdot\mathbf{B}_0)]=-\lambda\left(uB_0^2,0,0\right)
\end{equation}

Scaling the film thickness with the characteristic height of the flow
($h = h_0\tilde{h}$, $\delta = h_0 \tilde{\delta}$), the coordinates
by the characteristic length of the domain
$(x,z)=L(\tilde{x}, \epsilon \tilde{z})$ and the
velocity $(u, w) = (\nu/h_0\tilde{u}$, $\epsilon \nu/h_0 \tilde{w})$,
$U= (\nu/h_0)\tilde{U}$, the time $t=(h_0^2/\epsilon \nu)\tilde{t}$, the pressure
with $p = p_{a} + (\rho \nu^2/\epsilon h_0^2)\tilde{p}$ and the temperature
$T = T_{a} + \left(T_{s_0} - T_{a}\right)\tilde{T}$,
 where $\epsilon = h_0/L$ is the aspect ratio, $\nu = \mu/\rho$ is the
 kinematic viscosity of the fluid, $T_{s_0}$ is the temperature of the sheet at
the origin, and using the constitutive
relation \eqref{eq:cst} with Eqs.~\eqref{eq:r-e-tensor}, \eqref{eq:cstres}, and
 the external force \eqref{eq:efd},
the non-dimensional form of the governing
\eqref{eq:contd}-\eqref{eq:ic}, after dropping the tilde symbol in
explicit form are:
\begin{gather}\label{eq:contnd}
 u_x+w_z=0 \\
\label{eq:momndx}
\begin{aligned}
&\epsilon \left(u_t+uu_x+wu_z\right)  \\
&= -p_x+\epsilon^2 u_{xx}+u_{zz}+
 \mathrm{K}\bigl(\epsilon^3 u_{xxt}+ \epsilon u_{zzt}  +
 \epsilon^3uu_{xxx} -  \epsilon u w_{zzz} \\
&\quad + \epsilon^3 u_x u_{xx} -\epsilon u_xu_{zz} +\epsilon^3wu_{xxz}
 +\epsilon wu_{zzz} -  4\epsilon^3w_xw_{zz}-2\epsilon^5w_xw_{xx} \\
&\quad +  \epsilon^3 u_z w_{xx} -\epsilon u_zw_{zz} +2\epsilon w_zu_{zz}\bigr)
-\epsilon M^2 u\,,
\end{aligned} \\
\label{eq:momndz}
\begin{aligned}
&\epsilon^3 \left(w_t+uw_x+ww_z\right) \\
&= -p_z+\epsilon^4w_{xx}+\epsilon^2w_{zz}+
 \mathrm{K}\epsilon\bigl(\epsilon^4w_{xxt}+\epsilon^2w_{zzt}
 + \epsilon^2 \,u \,w_{xzz} \\
&\quad +  \epsilon^4\, u\, w_{xxx} +  2\epsilon^4\,u_x\,w_{xx}
 +\epsilon^2\,w_x\,u_{zz}-  \epsilon^4 \,w_x\, u_{xx} + \epsilon^2\, w\, w_{zzz} \\
&\quad -  \epsilon^4wu_{xxx}+\epsilon^2w_zw_{zz}-\epsilon^4w_zw_{xx}
-4\epsilon^2u_zu_{xx} -2u_zu_{zz}\bigr)-\epsilon \mathrm{Fr},
 \end{aligned} \\
\label{eq:tempnd}
\epsilon Pr \left(T_t+uT_x+wT_z\right) =\epsilon^2T_{xx}+T_{zz}, \\
\label{eq:bc2nd}
\begin{aligned}
&\epsilon S h_{xx} \left(1- M_w Ca T\right) (\epsilon^2 h_x^2 + 1)^{-1/2}\\
& =  - (\epsilon^2 h_x^2 + 1)p + 2\epsilon^2 \bigl(\epsilon^2 h_x^2 u_x
 - \epsilon^2 h_x w_x - h_xu_z+w_z\bigr) \\
&\quad + K\epsilon^3\bigl(2\epsilon^2u_{tx}h_{x}^2 -2\epsilon^2h_xw_{tx}-
 2h_xu_{tz}+2w_{tz}\bigr)\\
&\quad + K\epsilon\Bigl(\epsilon^2 h_x^2\bigl(2\epsilon^2uu_{xx}
 + 2\epsilon^2\,wu_{xz} +  u_z^2-\epsilon^4\,w_x^2\bigr)
 -  2\epsilon^2\,h_x\bigl(\epsilon^2\,uw_{xx} +uu_{xz} \\
&\quad  + \epsilon^2 w w_{xz}+ w u_{zz} +  \epsilon^2\,u_xw_x
 - \epsilon^2w_xw_z + u_zw_z-u_xu_z\bigr) \\
&\quad + 2\epsilon^2uw_{xz}+2\epsilon^2ww_{zz}+\epsilon^4 w_x^2- u_z^2\Bigr),
\end{aligned} \\
\label{eq:bc3nd}
\begin{aligned}
&-\epsilon M_w\left(T_x+h_xT_z\right)\left(1+\epsilon^2h_x\right)^{1/2} \\
&=(\epsilon^2w_x+u_z)(1-\epsilon^2h_x^2)
  + 2\epsilon^2 h_x(w_z-u_x) + K\bigl((1-\epsilon^2h_x^2)(\epsilon^3w_{tx}
 +\epsilon u_{tz}) \\
&\quad  +2\epsilon^3h_x(w_{tz}-u_{tx})\bigr)
 +K\Bigl((1-\epsilon^2h_x^2)  \bigl(\epsilon^3uw_{xx}+ \epsilon uu_{xz}
 + \epsilon^3ww_{xz} \\
&\quad + \epsilon wu_{zz}+ \epsilon^3u_xw_x
 -\epsilon^3w_xw_z +\epsilon u_zw_z-\epsilon u_xu_z\bigr)
 +2\epsilon h_x(\epsilon^2uw_{xz} \\
&\quad +\epsilon^2ww_{zz} - \epsilon^2uu_{xx}-\epsilon^2wu_{xz}
 + \epsilon^4w_x^2-u_z^2)\Bigr),
\end{aligned}\\
\label{eq:bc4nd}
 h_t = w - u h_x, \\
\label{eq:bc5nd}
 T_z - \epsilon^2h_xT_x=-B\left(1+\epsilon^2h_x^2\right)^{1/2} T, \\
\label{eq:bcnd1}
 \text{and at }z=0: \quad u=U(x), \quad w=0, \; T=\theta(x),\quad
\text{where } \theta(x)=\frac{T_{s}-T_{a}}{T_{s_0}-T_{a}}
\end{gather}

The non-dimensional parameters are the second grade parameter
$K=\alpha_1 / \rho h_0^2$,  the Hartmann number
$M=\sqrt{\lambda B_0^2h_0L/\rho\nu}$,
the Froude number $Fr=gh_0^3/\nu^2$, and the
Prandtl number $Pr=\rho c_p\nu/k$.

The symbol $S$ stands for the non-dimensional surface tension parameter
defined as $S=\epsilon^2 \sigma_a h_0 / \rho \nu^2$.
Here $Ca=\epsilon^2/S$, $B=\alpha h_0/k$ and
$M_w=\sigma_a h_0\gamma\left(T_{S_0}-T_a\right)/\rho\nu^2$
are the Capillary number, Biot number,  and the effective Marangoni number,
respectively. Equations \eqref{eq:bc2nd} and \eqref{eq:bc3nd} are obtained
after using the expression for unit normal vector
 $\hat{\textbf{n}}=\big(-h_{x}/\sqrt{1+h_{x}^2}, \;1/\sqrt{1+h_{x}^2}\big)$ and
 the unit tangent vector
$\hat{\textbf{t}}=\big(1/\sqrt{1+h_{x}^2}, \;h_{x}/\sqrt{1+h_{x}^2}\big)$.

The initial conditions read:
\begin{equation}\label{eq:icnd}
u = 0, \quad w =0, \quad  h(x, 0) = 1 + \delta(x).
\end{equation}

\section{Long-wave approximation}\label{sec:sect3}

The derivation of the one-dimensional thin film equation is based on long wave theory.
 The stated asymptotic analysis uses the
techniques of \cite{dandapat2011thin} but extends it to incorporate the
complex rheological effects.

For the underlying velocity and pressure variables, regular power series
expansion in $\epsilon$ are set up.
To obtain the equation of the thin film from the above problem, we expand
the variables as follows:
\begin{equation}\label{eq:as1}
 (u,w,p)=(u_0,w_0,p_0)+\epsilon(u_1,w_1,p_1)+\epsilon^2(u_2,w_2,p_2)+\dots
\end{equation}
Then the temperature variable $T$ is expanded as
\begin{equation}\label{eq:as2}
T=T_0+\epsilon^{1-n}T_1+\dots
\end{equation}
where $0<n<1$.

As per the scaling, the Prandtl number is of
order $O(\epsilon^{-n})$, the Froude number $Fr$ and the second-grade parameter
$K$ are of  order $O(1)$, the surface tension
parameter $S$ is of order $O(\epsilon^2)$, and the capillary number $Ca$ is of order
$O(1)$.

Substituting the expansions \eqref{eq:as1} and \eqref{eq:as2} in the
dimensionless model equations \eqref{eq:contnd}-\eqref{eq:bcnd1},
the \emph{leading order} equations are
\begin{gather}\label{eq:zome1}
 \frac{\partial u_0}{\partial x}+\frac{\partial w_0}{\partial z} = 0, \\
\label{eq:zome2}
-\frac{\partial p_0}{\partial x}+\frac{\partial^2 u_0}{\partial z^2} = 0, \\
 -\frac{\partial p_0}{\partial z} = 0, \\
 -\frac{\partial^2 T_0}{\partial z^2} = 0
\end{gather}
The corresponding boundary conditions are
\begin{gather}\label{eq:zobc2}
 \text{at } z=h(x,t):\quad p_0 = 0, \quad \frac{\partial u_0}{\partial z}=0,\quad
 \frac{\partial T_0}{\partial z} = -B T_0;\\
\label{eq:zobc1}
 \text{at }  z=0:\quad u_0=U(x),\quad w_0 = 0, \quad T=\theta(x)
\end{gather}

Using the stream function $\psi_0 = zU(x)$ the solution of the leading
 order problem is
 \begin{flalign}\label{eq:zos1}
 u_0=U, \quad w_0 = - zU_x, \quad p_0=0, \quad
 T_0= \Big(\frac{1+B(h-z)}{1+Bh}\Big)\theta(x)\,.
\end{flalign}
For the next order problem, we first solve the temperature equation.
With the results of zeroth-order, the temperature equation at order
$O(\epsilon^{1-n})$ reads
\begin{equation}\label{eq:fote}
\frac{\partial T_0}{\partial t}+u_0\frac{\partial T_0}{\partial x} +
 w_0\frac{\partial T_0}{\partial z}=\frac{\partial^2 T_1}{\partial z^2}
\end{equation}
with boundary conditions
\begin{equation}\label{eq:fobcte}
 \text{at } z=0:\quad T_1=0, \quad\text{ and at }
 z=h(x,t):\quad \frac{\partial T_1}{\partial z} = -B T_1.
\end{equation}
The solution of \eqref{eq:fote} and \eqref{eq:fobcte} satisfies
\begin{equation}\label{eq:solT1}
\begin{aligned}
T_1&=\frac{U\theta_x}{2}\Bigl\{z^2-zh\Bigl(\frac{2+Bh}{1+Bh}\Bigr)\Bigr\} \\
&\quad  - \Big\{\frac{B^2\theta U_xh}{\left(1+Bh\right)^2}+
\frac{B\left(U\theta_x-U_x\theta\right)}{1+Bh}\Big\}
\Big\{\frac{z^3}{6}-
\frac{zh^2}{6}\Big(\frac{3+Bh}{1+Bh}\Big)\Big\}
\end{aligned}
\end{equation}
The results \eqref{eq:solT1} is similar to the one derived in \cite{santra2009}
for the Newtonian fluid.

Using the results of zeroth-order and of order $O(\epsilon^{1-n})$, the final
solution of the temperature field reads
\begin{equation}\label{temperature_final}
\begin{aligned}
T&= \Bigl(\frac{1+B(h-z)}{1+Bh}\Bigr)\theta(x)+
\epsilon Pr \Bigl[\frac{U\theta_x}{2}\Bigl\{z^2-zh\Bigl(\frac{2+Bh}{1+Bh}\Bigr)\Bigr\}
 \\
&\quad - \Big\{\frac{B^2\theta U_xh}{\left(1+Bh\right)^2}+\frac{B\left(U\theta_x-
U_x\theta\right)}{1+Bh}\Big\}\Big\{\frac{z^3}{6}
-\frac{zh^2}{6}\Big(\frac{3+Bh}{1+Bh}\Big)\Big\}\Big]
\end{aligned}
\end{equation}

To close the problem, for the derivation of the free surface equation, we need
to solve the first and second order continuity and momentum equations.
The \emph{first order} problem reads:
 \begin{gather}\label{eq:contfo}
 \frac{\partial u_1}{\partial x}+\frac{\partial w_1}{\partial z}=0, \\
\label{eq:momfox}
\begin{aligned}
&\frac{\partial u_0}{\partial t}+u_0\frac{\partial u_0}{\partial x} +
 w_0\frac{\partial u_0}{\partial z} \\
& = -\frac{\partial p_1}{\partial x}
	+\frac{\partial^2 u_1}{\partial z^2}
	+K\Bigl[\frac{\partial^3 u_0}{\partial t\partial z^2}
-u_0\frac{\partial^3 w_0}{\partial z^3}
	- \frac{\partial u_0}{\partial x}\frac{\partial^2 u_0}{\partial z^2} \\
&\quad +  w_0\frac{\partial^3 u_0}{\partial z^3}
 -\frac{\partial u_0}{\partial z}\frac{\partial^2 w_0}{\partial z^2}
+2\frac{\partial w_0}{\partial z}\frac{\partial^2 u_0}{\partial z^2}\Bigr]-M^2u_0,
\end{aligned} \\
\label{eq:momfoz}
 -\frac{\partial p_1}{\partial z}-2K\frac{\partial u_0}{\partial z}
\frac{\partial^2 u_0}{\partial z^2}-Fr=0
 \end{gather}
with boundary conditions at the free surface $z=h(x,t)$:
 \begin{gather}\label{eq:bcm2}
 -p_1- K\left(\frac{\partial u_0}{\partial z}\right)^2 = S h_{xx}, \\
\label{eq:bcm3}
\begin{aligned}
&\frac{\partial u_1}{\partial z}+K\Bigl[\frac{\partial^2 u_0}{\partial t\partial z}+
u_0\frac{\partial^2 u_0}{\partial x\partial z}
 +w_0\frac{\partial^2 u_0}{\partial z^2}
 +\frac{\partial u_0}{\partial z}\frac{\partial w_0}{\partial z}\\
& - \frac{\partial u_0}{\partial x}\frac{\partial u_0}{\partial z} 
 - 2h_x\left(\frac{\partial u_0}{\partial z}\right)^2 \Bigr] \\
&= -M_w\Bigl[\frac{\partial T_0}{\partial x}+h_x\frac{\partial T_0}{\partial z}\Bigr]
 \end{aligned}
\end{gather}
and at the sheet $z=0$,
\begin{equation}\label{eq:bcm1}
 u_1=0,\quad  w_1=0
\end{equation}
Using the solution of zeroth order problem, the first order solutions for the
velocity and pressure can be found out as follows:
\begin{gather}\label{eq:fosm1}
 u_1=\left(UU_x+Frh_x-Sh_{xxx}+M^2U\right)\Big(\frac{z^2}{2}-hz\Big)+z\phi(x,t), \\
\label{eq:fosm2}
\begin{aligned}
 w_1
&=-\frac{z^3}{6}\left(UU_x+Frh_x-Sh_{xxx}+M^2U\right)_x \\
&\quad + \frac{z^2}{2}\frac{\partial}{\partial x}\bigl[h\left(UU_x+Frh_x-Sh_{xxx}
 +M^2U\right)-\phi(x,t)\bigr]
\end{aligned} \\
\label{eq:p1}
 p_1=Fr(h-z)-Sh_{xx}
\end{gather}
We introduce for brevity the notation
 $\phi(x,t)=-M_w\left(\frac{\theta}{1+Bh}\right)_x$ and
$f(x,t)= UU_x+Frh_x-Sh_{xxx}+M^2U$.

The model equations for \emph{second order} problem for velocity and pressure are:
\begin{gather}\label{eq:conso}
 \frac{\partial u_2}{\partial x}+\frac{\partial w_2}{\partial z}=0, \\
\label{eq:momsox}
\begin{aligned}
&\frac{\partial u_1}{\partial t}+u_0\frac{\partial u_1}{\partial x} +
 u_1\frac{\partial u_0}{\partial x}+w_0\frac{\partial u_1}{\partial z}
 +w_1\frac{\partial u_0}{\partial z} \\
&=-\frac{\partial p_2}{\partial x}+\frac{\partial^2 u_0}{\partial x^2}
 +\frac{\partial^2 u_2}{\partial z^2}\\
&\quad +  K\Bigl[\frac{\partial^3 u_1}{\partial t\partial z^2} -
 u_0\frac{\partial^3 w_1}{\partial z^3}- u_1\frac{\partial^3 w_0}{\partial z^3}
-\frac{\partial u_0}{\partial x}\frac{\partial^2 u_1}{\partial z^2} \\
&\quad -\frac{\partial u_1}{\partial x}\frac{\partial^2 u_0}{\partial z^2}
+w_0\frac{\partial^3 u_1}{\partial z^3}
+w_1\frac{\partial^3 u_0}{\partial z^3}
-\frac{\partial u_0}{\partial z}\frac{\partial^2 w_1}{\partial z^2}
-\frac{\partial u_1}{\partial z}\frac{\partial^2 w_0}{\partial z^2} \\
&\quad + 2\frac{\partial w_0}{\partial z}\frac{\partial^2 u_1}{\partial z^2}
 + 2\frac{\partial w_1}{\partial z}\frac{\partial^2 u_0}{\partial z^2}\Bigr]
 -M^2u_1\,,
\end{aligned}\\
\label{eq:momsoz}
 -\frac{\partial p_2}{\partial z}+\frac{\partial^2 w_0}{\partial z^2} +
 K\Big[-2\frac{\partial u_0}{\partial z}\frac{\partial^2 u_1}{\partial z^2}
-2\frac{\partial u_1}{\partial z}\frac{\partial^2 u_0}{\partial z^2}\Big] =0
\end{gather}
with boundary conditions at $z=h(x,t)$
\begin{gather}\label{eq:bcso2}
-p_2+2\Big(-h_x\frac{\partial u_0}{\partial z}+\frac{\partial w_0}{\partial z}\Big)
 +K\Big(-2\frac{\partial u_0}{\partial z}\frac{\partial u_1}{\partial z}\Big)=0,\\
\label{eq:bcso3}
\begin{aligned}
&\frac{\partial w_0}{\partial x}+\frac{\partial u_2}{\partial z}-
h_x^2 \frac{\partial u_0}{\partial z}+2h_x\Big(\frac{\partial w_0}{\partial z}
 -\frac{\partial u_0}{\partial x}\Big) \\
&  +K\Bigl[\frac{\partial^2 u_1}{\partial t\partial z}+
u_0\frac{\partial^2 u_1}{\partial x\partial z}
+ u_1\frac{\partial^2 u_0}{\partial x\partial z}
+  w_0\frac{\partial^2 u_1}{\partial z^2}
 +w_1\frac{\partial^2 u_0}{\partial z^2}
 +\frac{\partial u_0}{\partial z}\frac{\partial w_1}{\partial z} \\
& +\frac{\partial u_1}{\partial z}\frac{\partial w_0}{\partial z}
 - \frac{\partial u_0}{\partial x}\frac{\partial u_1}{\partial z}
 - \frac{\partial u_1}{\partial x}\frac{\partial u_0}{\partial z}
 + 2h_x\Big(-2\frac{\partial u_0}{\partial z}
 \frac{\partial u_1}{\partial z}\Big)\Bigr]=0,
\end{aligned}
 \end{gather}
and at $z=0$, i.e.
\begin{equation}\label{eq:bcso1}
 u_2=0, \quad w_2=0
\end{equation}
 Using the solution of the zeroth-order and first-order problem,
the second order problem  satisfies
\begin{gather}\label{eq:solso1}
\begin{aligned}
 u_2&=\Big(\frac{z^4}{24}-\frac{h^3z}{6}\Big)g_1(x,t)
+\Big(\frac{z^3}{6}- \frac{h^2z}{2}\Big)g_2(x,t)
 + \Big(\frac{z^2}{2}-hz\Big)g_3(x,t)  \\
&\quad  +  z\Bigl[U_{xx}h+4U_xh_x+Kf\left(h_t+Uh_x+U_xh\right) \\
&\quad +K\left(-\phi_t-U\phi_x+2 U_x\phi\right)\Bigr],
\end{aligned}\\
\label{eq:solso2}
\begin{aligned}
 w_2 &=\Big(-\frac{z^5}{120}+\frac{h^3z^2}{12}\Big)\frac{\partial g_1}
 {\partial x} + \frac{z^2h^2h_x}{4}g_1+\Big(-\frac{z^4}{24}
 +\frac{h^2z^2}{4}\Big)\frac{\partial g_2}{\partial x}\\
&\quad  +\frac{z^2hh_x}{2}g_2 
 + \Big(-\frac{z^3}{6}+\frac{hz^2}{2}\Big)\frac{\partial g_3}{\partial x}
 +\frac{z^2}{2}h_xg_3
 -\frac{z^2}{2}\frac{\partial}{\partial x}\Bigl[U_{xx}h \\
&\quad +4U_xh_x  + Kf\bigl(h_t + Uh_x +  U_xh\bigr) 
 +  K\left(-\phi_t-U\phi_x+2 U_x\phi\right)\Bigr],
\end{aligned} \\
\label{eq:solso3}
 p_2=-2U_x,
\end{gather}
where $g_1(x,t)=f_t+Uf_x-U_xf+M^2f$,
$g_2(x,t)=-(h_t+Uh_x+M^2h)f-hf_t-Uhf_x+\phi_t+U\phi_x+ M^2\phi$, and
 $g_3(x,t)=-3U_{xx}-K\left(f_t+Uf_x-3U_xf\right)$.

\section{Thin film equation}\label{sec:sect4}


Using the kinematic boundary condition, the free surface
evolution equation can be obtained as follows:
\begin{equation} \label{eq:tf}
 h_t+\frac{\partial}{\partial x}F(h)=0,
\end{equation}
where
\begin{align*}
F(h)&=Uh+\epsilon\bigl(-\frac{h^3f}{3}+\frac{h^2\phi}{2}\bigr)
 +\epsilon^2 \Bigl[\frac{-3h^5}{40}(f_t+Uf_x-U_xf+M^2f) \\
&\quad -\frac{5h^4}{24}\bigl(-fh_t - f_th-Uhf_x-Uh_xf
 +\phi_t+U\phi_x+ M^2\phi-M^2fh\bigr) \\
&\quad -\frac{h^3}{3}\Bigl(-3U_{xx}-Kf_t-KUf_x   + 3KU_xf\Bigr)
 +\frac{h^2}{2}\big(U_{xx}h+4U_xh_x \\
&\quad +Kf\left(h_t+Uh_x+U_xh\right) 
 + K\left(-\phi_t-U\phi_x+2U_x\phi\right)\big)\Bigr]
\end{align*}
with $\phi(x,t)=-M_w(\frac{\theta}{1+Bh})_x$ and
$f(x,t)=UU_x+Frh_x-Sh_{xxx}+M^2U$.

The closure of the thin film equation requires the boundary condition for the
film height. At the origin, we apply the symmetry boundary condition which
imposes that the gradient of the field
variable $h$ and its higher derivatives must vanish i.e.
\begin{equation}\label{eq:tfbc}
 h_x=0,\quad h_{xxx}=0, \quad h_{xxxx}=0.
\end{equation}
At the other end of the domain, we assume that the same sheet
stretching rate continues beyond the computed domain. We
also assume that the gradient of the free surface extends out
of the computational domain. These boundary conditions are
consistent with those mentioned in \cite{dandapat2011thin, singh}.
Finally, the model equation is supported by the initial condition
$h(x,0)= 1 + \delta(x)$.

\section{Numerical solution}\label{sec:sect5}

The transient film thickness \eqref{eq:tf} with supportive boundary and initial
conditions are solved using finite volume method.
 In this paper, we mostly follow the finite volume
technique described in Sellier et al.~\cite{sellier, panda}
 on a uniform grid system with implicit flux discretization.

\begin{figure}[htb]
 \begin{center}
 \includegraphics[width=0.7\textwidth]{fig2} % fvg  ppsfvg}
\end{center}
 \caption{Typical grid used for the finite volume discretization}
 \label{fig:fvd}
\end{figure}

 As shown in Figure \eqref{fig:fvd}, we discretized the fluid domain
using a uniform grid system.
 The unknown variable $h$ and the known fields $U$ and $\theta$ are located
at the cell centers.
 The flow domain $[0, L]$ is discretized into N equal size grid cells
of size $\Delta x_i = L/N$, and define $x_i= \Delta x_i/2 + (i-1)\Delta x_i$,
 $i= 1, 2, \dots, N$, so that $x_i$ is the center of the cell.
 The numerical solution is evaluated at
the discrete time levels $t^n$, $n=0,1,2,\dots$ with time step
$\Delta t^{n+1}=t^{n+1}-t^n$.
 From the given cell averaged solution $h_i^{n}$ over the cell $[x_i,x_{i+1}]$,
 the solution at the next  time step $t^{n+1}$ is obtained by
integrating \eqref{eq:tf}  over the space and time intervals
$[x_i,x_{i+1}]\times[t^n,t^{n+1}]$ which gives the discretized equation:
 \begin{equation}\label{eq:dtf}
 (h_i^{n+1}-h_i^n)\Delta x_i+(F_{i+1/2}^{n+1}-F_{i-1/2}^{n+1})\Delta t^{n+1}=0
 \end{equation}
for nodes $i=1,2,\dots,N$, where the discrete flux function is
\begin{equation}
 F_{i+1/2}^{n+1}=F(x_{i+1/2},t^{n+1}).
\end{equation}

 For the internal nodes, the face values are evaluated using
linear interpolation from nodal values and gradients using central differences
 such as
 \begin{gather*}
 h(x_{i+1/2},t^{n+1})=\frac{1}{2}\left( h_{i+1}^{n+1}+h_i^{n+1}\right),\\
 h_x(x_{i+1/2},t^{n+1})=\frac{1}{\Delta x_i} \left(h_{i+1}^{n+1}-h_i^{n+1}\right)
 \end{gather*}
Similarly, for other terms, we can derive the expressions like for the second,
third and fourth order derivative terms
 \begin{gather*}
 h_{xx}(x_{i+1/2},t^{n+1})
=\frac{1}{\Delta x_i^2} \left(h_{i+1}^{n+1}-2h_i^{n+1}+h_{i-1}^{n+1} \right), \\
 h_{xxx}(x_{i+1/2},t^{n+1})
=\frac{1}{\Delta x_i^3} \left(h_{i+2}^{n+1}-3h_{i+1}^{n+1}
+3h_i^{n+1}-h_{i-1}^{n+1} \right), \\
 h_{xxxx}(x_{i+1/2},t^{n+1})
=\frac{1}{\Delta x_i^4} \left(h_{i+3}^{n+1}-4h_{i+2}^{n+1}+6h_{i+1}^{n+1}
-4h_i^{n+1}+h_{i-1}^{n+1}\right )
 \end{gather*}

The boundary nodes at the origin are discretized using boundary conditions
\eqref{eq:tfbc} and at the other end of the domain; we
assume that the same sheet temperature profiles and the sheet
stretching rate continues beyond the computed domain. We
also assume that the gradient of the free surface extends out
of the computational domain.

Equation \eqref{eq:dtf} describes an implicit time discretization
scheme. Since the governing equation is nonlinear, a system
of nonlinear algebraic equations needs to be solved at each
time step. The MATLAB \text{\it{fsolve}} routine is used for this purpose.
This routine uses a non-linear least-squares algorithm to
solve a system of non-linear equations.
 A good initial starting guess is required to solve the nonlinear
equations. A reasonable initial guess for the free surface is
chosen to be unity throughout the discrete domain at the first
time step. Therefore, the disturbance function to the initial film
thickness $\delta(x)$ is considered to be zero throughout the discussion.
The solution from the previous time step can be
used otherwise. Convergence is usually achieved in less than
ten iterations, and the convergence criterion is that the norm
of the residuals should be less than $10^{-7}$.

\section{Results and discussion}\label{sec:sect6}

In this section some simulation results of free surface profiles and the
temperature distribution at the free surface are presented at different
times and for the different flow parameters.
In Figure \eqref{fig:hMv} the free surface profiles are given at different
times with variation of the second-grade non-Newtonian parameter.
The left panel results are for the fixed Hartmann number $M=2$ and the right
panel are for $M=4$. The other parameters used for the simulations are mentioned
in the figure caption. By comparing the left and right panels of
Figure \eqref{fig:hMv} we can see that with increasing
Hartmann number $M$ the fluid thinning behavior slows down due to the effect
of the Lorentz force.
It is further clear that the free surface height is decreasing with increasing
time. Again with increasing second-grade parameter the film thinning behavior
 decreases. This is due to the fact that for the higher value of $K$, the
tensile stress between fluid layer is higher and consequently higher resistance
to motion which decreases the rate of film thinning.

\begin{figure}[htb]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig3a} % h_t_v_k_v_M2_f
\includegraphics[width=0.7\textwidth]{fig3b} % h_t_v_k_v_M4_f
 \end{center}
 \caption{Film thickness profile for $S=2$, $Fr=2$, $B=1$, $M_{w}=2$,
 $\epsilon=0.1$ with
 $U=0.1x$, and $\theta(x)= e^{-x^2/33}$, left: Hartmann number $M=2$,
 right: Hartmann number $M=4$}
 \label{fig:hMv}
\end{figure}

The effect of Marangoni number $M_{w}$ on the free surface profile is described
in the left panel of Figure \eqref{fig:hMwv}. The sinusoidal
temperature profile $\theta(x)=0.6+0.5 \sin(2\pi x/10)$
(right panel of Figure \eqref{fig:hMwv}) is imposed
at the sheet and the linear stretching velocity $U=0.1x$ is considered.
The free surface profiles at two different times are
given for different values of $M_{w}$ with the variation of the second-grade
parameter $K$. It is clear from Figure \eqref{fig:hMwv} that the thermocapillary
deforms the free surface and this deformation is advected downstream
by stretching velocity. The Marangoni number that characterizes the relation
between the temperature dependent surface tension and viscous forces.
 The high Marangoni number means the lower the viscous forces and thereby
the fluid moves faster and the thinning rate increases.
It can also be observed that the thinning rate is faster
near the origin for lower value of viscoelastic parameter.

\begin{figure}[htb]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig4a} % h_t_Mw_v
\includegraphics[width=0.7\textwidth]{fig4b} % theta_sin
\end{center}
 \caption{Film thickness profile for $S=2$, $Fr=2$, $B=1$, $M=2$,
$\epsilon=0.1$ with  $U=0.1x$, left: Film thickness profile ($K=0$ (blue curve),
$K=5$ (red curve), $K=10$ (black curve)),
 right: Corresponding temperature profile $\theta(x)=0.6+0.5 \sin(2\pi x/10)$}
 \label{fig:hMwv}
\end{figure}

\begin{figure}[htb]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig5a} % h_U_var
\includegraphics[width=0.7\textwidth]{fig5b} % U_var
\end{center}
 \caption{Results for the different stretching velocity profile with
 temperature profile  $\theta(x) = 1-e^{-x^2/33}$ with
 $S=2$, $Fr=2$, $B=1$, $M=2$, $M_{w}=8$, $\epsilon=0.1$ at time $t=10$
 (non-dimensional), left:
 Film thickness profiles ($K=0$ (blue curve), $K=5$ (red curve),
 $K=10$ (black curve)),
 right: Corresponding stretching velocity profiles}
 \label{fig:hthin}
\end{figure}

The effect of stretching velocity on film height is discussed next.
The temperature profile $\theta(x) = 1-e^{-x^2/33}$ was imposed on the
sheet and the following three different stretching velocity distribution
is considered:
\begin{equation} \label{eq:Uv}
\begin{gathered}
U(x) = 0.6(0.1x + 0.01x^2),\quad \text{(parabolic concave)}, \\
U(x) = 0.1x, \quad \text{(linear)}, \\
U(x) = 0.75(1-(0.1x-1)^2), \quad \text{(parabolic convex)}.
\end{gathered}
\end{equation}
The free surface profiles are given at a fixed time $t=10$ (non-dimensional).
The results confirm that the parabolic concave stretching velocity profile
makes the film thinner away from the origin. It is further observed that
influence of second grade parameter decelerate the film thinning behaviour
of the fluid.

To explore the effect of stretching rate on the free surface profile
 with variation of second-grade parameter, the finite
volume scheme was run for the two different linear stretching
velocity, i.e., $U=\eta x$ where $\eta \in \{0.05, 0.1\}$, and
varying the values of $K$ from $K=0$ to $K=10$
with an increment of five.
The free surface profile is plotted in Figure \eqref{fig:hetaK} at a
fixed time $t=10$ (non-dimensional)
for different stretching speed with different values of $K$.
For $U = \eta x$ where $\eta = 0.1$, the rapid stretching of the film will
result in the rapid build-up of stresses. Clearly, the faster the stretching,
the faster the thinning of the fluid.
The larger the viscoelastic parameter $K$, the larger the departure from the
purely viscous case ($K=0$)  because the additional stresses the fluid needs
to overcome the flow.

\begin{figure}[htb]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig6} % h_etax_f
\end{center}
 \caption{Effect of stretching rate and viscoelastic parameter $K$ on free surface
 profile with
 $S=2$, $Fr=2$, $B=1$, $M=2$, $M_{w}=8$, $\epsilon=0.1$ at time
$t=10$ (non-dimensional) with
 $\theta(x) = 1-e^{-x^2/33}$}
 \label{fig:hetaK}
\end{figure}

\begin{figure}[htb]
\begin{center}
\includegraphics[width=0.7\textwidth]{fig7a} % temp
\includegraphics[width=0.7\textwidth]{fig7b} % h_cor_for_T
\end{center}
 \caption{Variation of temperature profile at the interface ($z=h$ for different
 Biot number $B$ for $S=2$, $Fr=2$, $K=1$, $M=2$, $M_{w}=2$, $\epsilon=0.1$,
$Pr=1$ with  $U(x)= 0.6(0.1x+0.01x^2)$ and $\theta(x)=0.6 + 0.5 \sin(2 \pi x/10)$,
left: temperature profiles
 at the interface, right: corresponding free surface profiles}
 \label{fig:tempB}
\end{figure}

The effect of Biot number on the temperature profile at the free surface of the
 thin film is demonstrated in Figure \eqref{fig:tempB}. The temperature profile
at the fluid-gas interface is given at two different times with variation of
Biot number ($B$).
 The Biot number, that compares the relative magnitudes of resistances to
internal conduction and surface convection. In a high Biot number flow,
the surface convection dominates the internal conduction and consequently
the rate of heat transfer to the fluid through the sheet is higher.
Moreover, temperature increases with increase of time as
film thins faster (right panel of Figure \eqref{fig:tempB}).
We have observed that the Biot number has
no significant effect on the thinning behavior of the fluid.

\subsection*{Conclusion} %\label{sec:sect7}
A model for the dynamics of thin film over the heated steady stretching sheet
under the uniform transverse magnetic field is derived based on the
long-wave theory. The derivation of the numerical scheme for the unsteady
free surface model by finite volume method and their simulations for the
wide range of parameters provides a deeper understanding of the film thinning
process of the non-Newtonian second-grade fluid. The numerical solution reveals
the dependency of the film thinning behavior on the stretching rate and the
second-grade fluid parameter.
Moreover, the film thinning rate is influenced by the thermocapillary forces
and magnetic parameter.
It is hoped that the desired film thickness can be achieved by controlling
the flow parameters. Future work will explore the derivation of free surface
model for the dynamics of film thickness in a similar flow geometry
for a third-grade non-Newtonian fluid.
		
\begin{thebibliography}{00}

\bibitem{andersson1996flow} H.~I.~Andersson, J.~B.~Aarseth, N.~Braud,  B.~S.~Dandapat;
 Flow of a power-law fluid film on an  unsteady stretching surface, \emph{J. of
 Non-Newtonian Fluid Mechanics}, 62(1) (1996), 1--8.

\bibitem{Andersson2000} H.~I.~Andersson, J.~B.~ Aarseth, B.~S.~Dandapat;
 Heat transfer in a liquid film on an unsteady stretching
surface, \emph{Int.~J.~Heat Mass Transfer}, 43 (2000), 69--74.

\bibitem{Aziz2010} R.~C.~Aziz, I.~Hashim;
Liquid film on unsteady stretching sheet with general surface temperature and
viscous dissipation, \emph{Chin.~Phys.~Lett.}, 27 (2010), 110202.

\bibitem{Liu2008} I.-Chung~ Liu, H.~I.~Andersson;
 Heat transfer in a liquid film on an unsteady stretching sheet,
 \emph{Int.~J.~Therm. Sci.}, 47 (2008), 766--772.

\bibitem{dandapatZAMP2006} B.~S.~Dandapat, A.~Kitamura, B.~Santra;
Transient film profile of thin liquid film flow on a stretching surface,
\emph{Zeitschrift fur angewandte Mathematik
und Physik ZAMP}, 57(4) (2006), :623--635.

\bibitem{dandapatPoF2006} B.~S.~Dandapat, S.~Maity;
 Flow of a thin liquid film on an unsteady stretching sheet,
\emph{Physics of Fluids}, 18(10) (2006), Article ID 102101, 7 pages.

\bibitem{dandapat} B.~S.~Dandapat, B.~Santra, H.~I.~Andersson;
 Thermocapillarity in a thin liquid film on an unsteady stretching
surface, \emph{Int.~J.~Heat and Mass Transfer}, 46 (2003), 3009--3015.

\bibitem{vajravelu} B.~S.~Dandapat, B.~Santra, K.~Vajravelu;
The effects of variable fluid properties and
thermocapillarity on the flow of a thin film on an unsteady stretching
sheet, \emph{Int.~Jr.~Heat and Mass Transfer}, 50 (2007), 991--996.

\bibitem{dandapat2011thin} B.~S.~Dandapat, S.~K.~Singh;
 Thin film flow  over a heated nonlinear stretching sheet in
presence of uniform transverse magnetic field,
\emph{Int. Commu. in Heat and Mass Transfer}, 38(3) (2011), 324--328.

\bibitem{dunn} J.~E.~Dunn, R.~L.~Fosdick;
 Thermodynamics, stability and boundedness of fluids of
complexity 2 and fluids of second grade,
\emph{Arch Ration Mech Anal.}, 56(3) (1974), 191--252.

\bibitem{dunnetal} J.~E.~Dunn,  K.~R.~Rajagopal;
Fluids of differential type critical review and thermodynamics
analysis, \emph{Int. Jr. Eng Sci.}, 33:689--747, 1995.

\bibitem{fosdicketal} R.~L.~Fosdick, K.~R.~Rajagopal;
 Anomalous features in the model of second order fluids.
 \emph{Arch Ration Mech Anal.},  70:145--152, 1979.

\bibitem{hayat2008influence} T.~Hayat, S.~Saif, Z.~Abbas;
 The influence of heat transfer in an MHD second grade fluid film over
an unsteady stretching sheet, \emph{Physics Letters A}, 372(30) (2008),5037--5045.

\bibitem{hayat2007analytic} T.~Hayat, M.~Sajid;
 Analytic solution for axisymmetric flow and heat transfer of a second grade
fluid past a stretching sheet, \emph{Int.~ J.~of Heat and Mass Transfer},
 50(1) (2007), 75--84.

\bibitem{levy} R.~Levy, D.~B.~Hill, M.~G.~Forest, J.~B.~Grotberg;
Pulmonary fluid flow challenges for experimental and
mathematical modeling. \emph{Integrative and Comparative Biology},
54(6) (2014), 985--1000.

\bibitem{Khan2011} Y.~Khan, Q.~Wu, N.~Faraz, A.~Yildirim;
 The effects of variable viscosity and thermal conductivity on a thin film over
a shrinking/stretching sheet, \emph{Comp.~Math.~Appl.}, 61 (2011), 3391--3399.



\bibitem{ali} F.~Md. Ali, R.~Nazar, N.~Md.Arifin, I.~Pop;
MHD boundary layer flow and heat transfer over a stretching sheet
with induced magnetic field, \emph{Heat and mass Transfer},
47(2) (2011), 155--162.

\bibitem{Noor2010} N.~F.~M.~Noor, I.~Hashim;
 Thermocapillary and magnetic field effects in a thin liquid film on an
unsteady stretching  surface, \emph{Int.~Jr.~Heat Mass Transfer},
53 (2010) 2044--2051.

\bibitem{panda} S.~Panda, M.~ Sellier, M.~C.~S.~Fernando, M.~K.~Abeyratne;
 Process parameter identification in thin film flows driven by a stretching surface,
{Int. Jr. of Engineering Mathematics}, Article ID 485431, 2014.

\bibitem{rivlin} R.~S.~Rivlin, J.~L.~Ericksen;
 Stress deformation relations for isotropic materials,
 \emph{J. Ration Mech. Anal.}, 4 (1955), 323--425.

\bibitem{sajid2009unsteady} M.~Sajid, I.~Ahmad, T.~Hayat, M.Ayub;
 Unsteady flow and heat transfer of a second grade fluid over a stretching sheet,
 \emph{Comm. in Nonlinear Science and Numerical
Simulation}, 14(1) (2009), 96--108.

\bibitem{sajid2008influence} M.~Sajid, T.~Hayat;
 Influence of thermal radiation on the boundary layer flow due to an exponentially
stretching sheet, \emph{Int. Comm. in Heat and Mass Transfer},
 35(3):347--356, 2008.

\bibitem{salleh} M.~Z.~Salleh, R.~Nazar, I.~Pop;
 Boundary layer flow and heat transfer over a stretching
sheet with Newtonian heating, \emph{J. of the Taiwan
Institute of Chemical Engineers}, 41(6) (2010), 651--655.

\bibitem{santra2009} B.~Santra, B.~S.~Dandapat;
 Unsteady thin-film flow over a heated stretching sheet,
\emph{Int. Jr. of Heat and Mass Transfer}, 52 (2009), 1965--1970.

\bibitem{sellier} M.~Sellier, S.~Panda;
Surface temperature reconstruction based on the thermocapillary effect,
\emph{The ANZIAM Journal}, 52:146--159, 2010.

\bibitem{singh} S.~K.~Singh, B.~S.~Dandapat;
 Thin film flow of Casson liquid over a non-linear stretching sheet in the
presence of a uniform transverse magnetic field, \emph{Canadian Journal of Physics}
93 (2015), 1067--1075.

\bibitem{Wang} B.~Wang;
 Liquid film on an unsteady stretching surface,
 \emph{Quart.~Appl.~Math.}, 48 (1990), 601--610.

\end{thebibliography}

\end{document}



